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ABSTRACT 


One of the primary advantages of laminated composites is the ability to orient the fibers of 
each layer to optimize its mechanical properties such as strength and convenient stiffness. 
Notwithstanding these excellent properties, laminates have a problem specific to composite 
materials made by lamination: This is the inter-laminar rupture, known under the technical 
name by the 'Delamination'. The attempt to produce functional gradient materials (FGM) was 
born from the mixing of two materials in such a way that the two materials preserve their 
physical, mechanical and thermal properties most effectively. The resultant material FGM 
shows a gradation across the depth of the structure, from a typically metallic material such as 
steel or aluminum on one side of the depth, to another material such as ceramic on the other 
side. Thus, the functional gradient or functionally graded material (FGM) is a new variety of 
composite materials having a gradual and continuous variation of the volume fractions of each 
constituents (in general, metal and ceramic) across the thickness, inducing changes as a result 
of the global thermo-mechanical properties of the structural element it constitutes. The change 
in the properties of the materials for our study takes according to a power law. This allows 
these composites the capability to significantly reduce the delamination effect and stress 
concentration encountered in conventional composites, but the most important 1s to have a 
new function that did not exist before. Therefore, the FGM material designed to overcome the 
inter-laminar rupture problem (delamination) on the one hand and on the other hand to 
remedying the problems caused by severe thermal environments. The analysis of structures 
in FGM, requires to establish a tools for modeling mechanical behavior increasingly 
sophisticated, including finite element calculation 1s essential for conception, designated and 
verification of these complex structures. In the present thesis, we develop a new transverse 
shear function based on a high order beam theory to analyze the bending behavior of FGM 
beams. The present theory takes into account the transverse shear effect on the one hand and 
on the other hand, it does not require the use of transverse shear correction factor. It gives a 
parabolic description of the shear stress across the thickness while respecting the condition of 
zero shear stress on the free edges of the FGM structure (beam). We made a comparison 
between the different high orders theories with that developed in the present work. This 
comparison highlighted the present shear function with those found in the literature. Then, we 
develop a new flat shell element based on the strain formulation for the case of an isotropic 


material. The element is the result of the superposition of a membrane element with a bending 


plate element. It is a plane triangular element with three nodes and six degrees of freedom for 
each node. Therefore, the rotation around the normal is considered. The comparison of our 
results with reference solutions shows the good performance and precision of the formulated 
shell element. In a second time, we enrich and enlarge the development of this shell element 
so that it can model FGM structures. For this purpose, we extend additional developments for 
the established element to model FGM structures. In this case, a membrane-bending coupling 
is necessarily to be taken into account. The proposed element provides robust and accurate 
linear analyzes of FGM shells relative to the displacement approach. Comparing our results 
with confrontational solutions shows that the developed element is competitive, robust, 


efficient and accurate. 


Key words: Material, composite, gradient, functional, beam, membrane, plate, shell, finite 


element, strain formulation, static condensation, drilling rotation, and delamination. 


RESUME 


L’un des avantages primordiaux des composites stratifiés est l’aptitude d’orienter les 
fibres de chaque pli afin d’optimiser ses propriétés mécaniques telles que la résistance et la 
rigidité utile. Nonobstant ces excellentes propriétés, les stratifiés présentent un problème 
propre aux matériaux composites réalisés par stratification : Il s’agit de la rupture inter- 
laminaire, connue sous l’appellation technique par le ‘Délaminage”. La tentative de produire 
des matériaux à gradient fonctionnel (FGM) est née du mélange de deux matériaux de telle 
manière que les deux matériaux conservent leurs propriétés physiques, mécaniques et 
thermiques le plus efficacement. Le matériau résultant FGM montre une gradation à travers la 
profondeur de la structure, d'un matériau typiquement métallique tel que l'acier ou l'aluminium 
à une face de la profondeur, à un autre matériau tel que la céramique sur l'autre face. Donc, le 
matériau à gradient fonctionnel ou fonctionnellement gradué (FGM) est une nouvelle gamme 
de matériaux composites ayant une variation graduelle et continue des fractions volumiques 
de chacun des constituants (en général, métal et céramique) à travers l’épaisseur, induisant des 
changements, en conséquence des propriétés thermomécaniques globales de l’élément 
structural qu’il constitue. Le changement des propriétés des matériaux est pris selon une loi de 
puissance pour notre étude. Cela permet à ces composites la possibilité de réduire 
significativement l’effet du délaminage et de la concentration des contraintes rencontrées chez 
les composites conventionnels, mais le plus important est d'avoir une nouvelle fonctionnalité 
qui n'existait pas avant. Donc, le matériau FGM a été conçus pour pallier au problème de la 
rupture inter-laminaire (le délaminage) d’une part et pour remédier aux problèmes engendrés 
par des environnements thermiques sévères d’autre part. L’analyse des structures en FGM 
nécessite de mettre en place des outils de modélisation du comportement mécanique de plus 
en plus sophistiqués, notamment, le calcul par la méthode des éléments finis est indispensable 
pour le dimensionnement et la vérification de ces structures complexes. Dans cette thèse, une 
nouvelle fonction de cisaillement transversal basé sur une théorie de poutre d’ordre élevé a été 
développée afin d’analyser le comportement en flexion des poutres FGM. La présente théorie 
prend en compte l’effet de cisaillement transverse d’une part et elle n'exige pas l’utilisation de 
facteur de correction de cisaillement transversal d’autre part. Elle donne une description 
parabolique de la contrainte de cisaillement à travers l’épaisseur tout en respectant la condition 
de contrainte de cisaillement nulle sur les bords libres de la structure (poutre) en FGM. Une 
comparaison entre les différentes théories d’ordre élevé avec celle développée dans le présent 


travail a été dressée. Cette comparaison a mis en valeur la présente fonction de cisaillement 


avec celles se trouvant dans la littérature. Puis, un nouvel élément de coque à facette plane 
basé sur la formulation en déformation a été développé pour le cas d’un matériau isotrope. 
L’élément est le résultat de la superposition d’un élément de membrane et d’un élément de 
plaque en flexion. C’est un élément triangulaire plan, à trois nœuds et à six degré de liberté 
par nœuds. Donc, la rotation autour de la normale est prise en considération. La comparaison 
de nos résultats avec des solutions de références montre la bonne performance et la précision 
de l’élément de coque formulé. Dans un deuxième temps le développement de cet élément de 
coque est enrichi et élargie pour qu’il puisse modéliser des structures en FGM. A cet effet des 
développements supplémentaires ont été étendu pour que l’élément puisse modéliser des 
structures en FGM. Dans ce cas un couplage membrane-flexion est forcément à prendre en 
compte. L'élément proposé fournit des analyses linéaires robustes et précises des coques en 
FGM par rapport à l'approche fondée sur le déplacement. La comparaison de nos résultats avec 
des solutions de comparaison montre que l’élément développé est compétitif, robuste, 


performant et précis. 


Mots clés : Matériau composite, gradient fonctionnel, poutre, membrane, plaque, coque, 


élément fini, formulation en déformation, condensation statique, drilling rotation, délaminage. 
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1. Généralités 

Les matériaux ont joué un rôle important dans la société humaine tout au long 
de l'histoire. L'humanité a toujours essayé de produire des matériaux plus résistants et 
très solides pour construire des structures durables servant d'abri pour l’être humain. 
Aux premiers temps de la civilisation (1500 av. J.C.), les Égyptiens et les 
Mésopotamiens mélangent de la paille et de la boue pour fabriquer des briques dans le 
but de construire des bâtiments assez résistants et plus durables. Par la suite, en 1800, 
le béton devient un composite largement utilisé. Il se génère en mélangeant le ciment, 
l'agrégat et l'eau. Plus tard, au début des années 1900, le plastique renforcé de fibres, 
un composite composé d'une matrice polymère renforcée de fibres, devient un matériau 
fondamental pour l'industrie aérospatiale, automobile, maritime et de la construction. 
Actuellement, la meilleure manière pour combiner deux matériaux, en conservant leurs 
propriétés, consiste à les assembler d'un pourcentage variable sur la section 
transversale. Par exemple, l’ancrage d'une surface en céramique sur le dessus et d'une 
surface en acier sur le dessous avec une zone de transition lisse, appelée matériau à 
gradation fonctionnelle (FGM). II est accepté comme un meilleur moyen de collage 
plutôt que de les déposer directement sans gradation. Par cette technique, les deux 
matériaux préservent leurs meilleures propriétés, comme la haute résistance, la rigidité 
élevée, la résistance à haute température, la faible densité de la céramique et la ténacité 
et la malléabilité de l'acier. Les FGM s’utilisent dans plusieurs industries et secteurs 
comme l'aérospatiale, le nucléaire, la défense, l'automobile, la communication, 
l'énergie, etc. En plus d'être produits artificiellement, les formes primitives des FGM 
existent dans la nature. L’os, la peau humaine, l’arbre de bambou peuvent être considéré 
comme des formes physiologiques d’FGM. La figure 1 illustre plusieurs exemples 
organiques et artificiels évoquant les FGM. 
Les domaines d'application des FGM se résument comme suit : 
L’Aérospatial : Les boucliers thermiques pour engins spatiaux, les tubes échangeurs de 
chaleur. 
Le Biomédical : L’os artificiel, la peau et les dents. 
La Communication : Les fibres optiques, les lentilles et les semi-conducteurs. 
Le Domaine du nucléaire : Les palettes de combustible, la paroi de plasma des réacteurs 


de fusion. 
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Le Secteur de l'énergie : Les générateurs thermoélectriques, les cellules solaires et les 


capteurs. 
L’ Automobile : Les systèmes de transmission de puissance et les systèmes de freinage. 
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Figure 1 Illustrations organiques et artificielles pour FGM 
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Le concept théorique de production des FGM a été proposé en 1984 au Japon 
(Jha et al., 2013). Le processus de fabrication et de conception des FGM devient le sujet 
capital d'intérêt au cours des trois dernières décennies. La fabrication des FGM se 
discute sous deux aspects, à savoir la gradation et la consolidation. La gradation forme 
la structure spatialement inhomogène. Les processus de gradation se catégorisent en 
processus constitutifs, d'homogénéisation et de ségrégation. Les processus constitutifs 
découlent d'une génération par étapes du corps gradué à partir de matériaux pionniers 
ou de poudres. Les progrès de la technologie d'automatisation au cours des dernières 
décennies ont assuré la viabilité technologique et économique des processus 
constitutifs. L'homogénéisation est la mutation de la transition subite entre deux 
matériaux en un gradient. Le processus de ségrégation commence avec un matériau 
macroscopiquement homogène, qui se transforme en un matériau gradué par le 
transport de matière provoqué par un champ externe (un champ magnétique ou 
gravitationnel). Les processus d'homogénéisation et de ségrégation produisent des 
gradients continus ; cependant de tels procédés ont des limites concernant la 
distribution des matériaux à produire. Le processus de consolidation, tels que le 
séchage, le frittage ou la solidification, vient généralement après le processus de 
gradation. Les critères de traitement doivent être déterminés de telle sorte que le 
gradient ne soit pas détruit ou modifié de façon non restreinte (Kieback et al. 2003). 
2. Enquête sur la modélisation par éléments finis des FGM 

Dans le but d’évoluer vers une production et une utilisation plus efficace des 
matériaux FGM dans les applications d'ingénierie, plusieurs techniques de conception, 
d'analyse et de calcul ont également été proposées. 

L’analyse des structures en FGM a connu un envol en utilisant les méthodes 
numériques surtout la méthode des éléments finis. 

Les poutres, les membranes, les plaque et les coques constituent des éléments de 
base dans les structures aérospatiales, marines et des construction civil et industriels, 
c’est pourquoi un intérêt remarquable est porté, pour bien appréhender leur 
comportement sous diverse sollicitations est une phase décisive dans l’analyse 
structurelle. 

La plupart des analyses des éléments structuraux en FGM reposent sur la théorie 
de l'élasticité. Il convient de mentionner que ces analyses proposent des solutions 
spécifiques aux problèmes de charge et aux conditions aux limites et elles se focalisent 


uniquement sur les comportements élastiques linéaires des matériaux. Une alternative 
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à une telle analyse serait le développement d'éléments finis : de poutre, de membrane, 
de plaque, de coque et d’éléments finis solides permettant l'analyse des problèmes de 
modélisation des structures en FGM. Il est à noter que même les logiciels disponibles 
comme NASTRAN, ANSYS ou ABAQUS ne fournissent pas automatiquement de 
telles capacités et nécessitent le développement de subroutines et de fonctions définies 
par l'utilisateur. Comme on peut le voir sur la figure 1 des exemples concernant les 
domaines d'utilisation des FGM (revêtement thermique, écrans thermiques, carters de 
fusées, etc.), les éléments de coque et de plaque sont plus adaptables et préférables pour 
la modélisation de telles structures. 

Une autre approche de modélisation par éléments finis par rapport à l'élément de 
coque ou de plaque est le développement d'éléments finis de poutre. Les théories des 
poutres décrivent la cinématique des sections d'un élément de poutre. Ces théories de 
poutres se regroupent sous les théories des poutres d'Euler-Bernoulli (EBBT Euler- 
Bernoulli Beam Theory), de Timoshenko (TBT : Timoshenko Beam Theory) et d'ordre 
supérieur (HOBT : High Order Beam Theory) dans la littérature. Le développement 
d'éléments finis de poutres par l'utilisation de diverses théories de poutres est 
principalement classé en deux, où le choix le plus largement adopté est l'utilisation de 
l'approche par déplacement et la méthode alternative est l'utilisation d'une approche de 
formulation en déformation. Dans le contexte des formulations basées sur le 
déplacement et des hypothèses de la théorie des poutres de Timoshenko pour des 
poutres prismatiques uniformes avec des matériaux homogènes, les éléments finis de 
poutre proposés par Friedman et Kosmatka (1993) et Reddy (1997) s'avèrent être un 
élément libre de tous blocage de cisaillement transversal et fournissent une matrice de 
raideur exacte sous une réponse élastique linéaire. 

En termes de développement de poutres en FGM, le travail le plus ancien et le 
plus cité est l'élément fini de poutre proposé par Chakraborty et al. (2003) dans lequel 
ils ont développé un nouvel élément de poutre pour étudier les interactions 
thermoélastiques des poutres FGM. La génération des fonctions de forme, de rigidité et 
de la matrice de masse pour une poutre FGM est expliquée dans cet article de recherche. 
La matrice de rigidité est dépourvue de verrouillage (blocage) de cisaillement 
transversal puisque le degré des fonctions de forme du déplacement transversal est plus 
élevé que la fonction de forme de la pente comme cela est discuté par Reddy (1997). II 


convient de signaler que l'étude de Chakraborty et al. (2003) est couramment cité dans 


Introduction générale 


la littérature et elle est considérée comme une formulation de référence pour la 
modélisation des poutres en FGM. 

Dans l'étude de Hemmatnezhad et al. (2013), l'analyse des vibrations libres de 
grande amplitude, imposant une grande déformation de type Von-Karman, de poutres 
avec des matériaux graduellement fonctionnels, est réalisée. La formulation d'éléments 
de poutre finie pour FGM, proposée par Chakraborty et al. (2003) est encore utilisé 
dans cette étude. Les résultats obtenus à partir des simulations numériques sont 
comparés aux solutions de la théorie d'élasticité disponibles dans la littérature. Outre le 
développement et l'utilisation d'éléments finis basé sur la formulation en déplacement, 
il est également possible de formuler des éléments finis basé sur la formulation en 
déformation ou la formulation mixte pour l'analyse des poutres en FGM. 

Les éléments de poutre restent limités dans leurs applications aux structures 

unidimensionnelles, ce qui impose la formulation d’éléments de coque et de plaque qui 
sont plus adéquates pour la modélisation de structures bidimensionnelles. 
A cet effet, Reddy et Chin en 1998 ont réalisé une étude sur le comportement mécanique 
et thermique des éléments de coque en FGM. Un élément fini de coque pour l'analyse 
des FGM est développé dans ce travail et des études paramétriques de l'effet de 
couplage thermomécanique sur les FGM sont également réalisées. Les paramètres sont 
sélectionnés en tant que rapports de distribution des matériaux et de leurs types. 

Les comportements thermomécaniques des plaques en FGM ont été analysés par 
Reddy en 2000. Il a utilisé des éléments finis basés sur une formulation de grande taille 
de type Von-Karman. Les résultats obtenus indiquent que la réponse de base des 
plaques en FGM constituées de métal et de céramique ne se situe pas forcément entre 
le métal total et l’état complètement en céramiques. La déflexion non-dimensionnelle 
atteint un minimum à un indice de fraction volumique dépendant des propriétés des 


constituants et de leurs rapports. 


3. Objectifs et portée de la thèse 

Le premier effort du présent travail consiste au développement d’une fonction de 
de cisaillement permettant de déterminer les contraintes de cisaillement transversal en 
chaque point de la section d’une poutre en FGM. La fonction développée vérifie les 
conditions de contraintes nulles sur les faces supérieure et inférieur de la poutre et celle 
de la contrainte maximale sur l’axe neutre avec un changement de pente aux points 


d’inflexions. Les résultats obtenus avec cette fonction sont très satisfaisants. 
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La deuxième étendue de la présente recherche concerne la construction d’un 
élément fini de coque mince basé sur la formulation en déformation permettant la 


modélisation et la discrétisation même des structures à doubles courbure. 


Le troisième effort de la présente recherche se concentre sur la généralisation de 
la formulation en déformation pour la modélisation des structures courbes en FGM par 
la construction d’un élément finis de coque dédié pour l’analyse des coques en FGM. 
Pour ce but des matrices de rigidité de membrane, de flexion et de couplage membrane- 
flexion sont développés afin de construire la matrice de rigidité de l’élément de coque 


en matériau FGM. L’élément a permis d’analyser quelques exemples test en FGM. 


Par ailleurs, un autre objectif de ce travail est d'étudier l'influence des différents 
paramètres tels que la géométrie, les conditions aux limites, l’épaisseur (mince où 
épaisse), le paramètre du matériau utilisé sur les différents types de comportements des 


structures en FGM. 


4, Organisation de la thèse 

Dans le premier chapitre nous faisons une présentation très large sur les matériaux 
composites : les composites classiques (stratifiés) en premier lieu, en montrant leurs 
limitations et les problèmes ouverts qui leurs sont inhérents et les matériaux composites 
en FGM en deuxième lieu, en révélant la nécessité de leur modélisation d’une manière 
adéquate. Donc, une présentation générale des matériaux composites, de leurs 
mécanismes de rupture sont ultérieurement abordés en focalisant sur le phénomène de 
délaminage, puis nous donnons une description des Matériaux à Gradient Fonctionnel 
‘MGF’, leur histoire d’évolution, leurs propriétés et leurs domaines d’utilisation. 

Au chapitre deux, nous présentons les relations de la cinématique et de l’équilibre 
nécessaires pour le développement de l’élément de poutre avec cisaillement transversal 
dédié pour les matériaux à gradient fonctionnel (FGM). Ainsi nous mentionnons 
quelques théories principales de poutres et les cinématiques de leur section. Ensuite, 
nous présentons une nouvelle fonction de cisaillement basé sur une théorie de 
déformation de cisaillement d’ordre élevé (HSDT : Higher Shear Deformation Theory) 
destinée à la modélisation des poutres avec cisaillement en matériau FGM. Un état 


d’analyse comparative sur la base d’une étude paramétrique des performances de 


E 
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l’élément développé en tenant en compte de la nouvelle fonction de cisaillement par 


rapport à ceux développés par d’autres auteurs est exposée à la fin du chapitre. 


Le troisième chapitre est réservé à la formulation et au développement de 
l’élément de coque mince triangulaire à facette plane à trois nœuds et à six degré de 
liberté par nœud basé sur le modèle en déformation pour le cas d’un matériau isotrope. 
La matrice de rigidité membranaire est celle de l’élément triangulaire de membrane à 
trois nœuds et à trois degré de liberté par nœud développé sous la formulation en 
déformation et basé sur les équations d’équilibre. La matrice de rigidité flexionnelle est 
celle de l’élément de plaque mince triangulaire à trois nœuds et à trois degré de liberté 
par nœuds développé sur la base du modèle en déformation. Une validation très riche 
est effectuée sur la base de tests de référence (Benchmark) pour les éléments finis de 


plaques et coques. 


Au quatrième chapitre nous généralisons l’application de l’élément formulé au 
chapitre trois pour qu’il soit capable de modéliser des coques en matériau FGM et cela 
par un développement complémentaire. En plus de l’adaptation des deux éléments de 
membrane et de flexion développés au chapitre trois pour les matériaux FGM, une 
matrice de rigidité de couplage membrane-flexion est indispensable dans le cas des 


FGM et cela suite à la non confusion du plan médian et du plan neutre de rigidité. 


Finalement, notre thèse s’achève par une conclusion générale mettant en valeur 


les résultats obtenus et dégageant des perspectives pour des travaux futurs. 
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Suite à la production originale d’un matériau composite moderne baptisé : matériau à 
gradient fonctionnel; ce chapitre est consacré à présenter l’aplomb du problème. Une 
introduction générale sur le matériau composite et le rôle des différentes matériaux le 
composant (matrice, renfort, charge, interface...) seront étudiés en mettant le point sur leurs 
avantages et leurs inconvénients. Les mécanismes de rupture des composites 
unidirectionnels seront ensuite abordés en insistant sur le phénomène de délaminage observé 
dans les composites multicouches. Ensuite, un historique détaillé du développement et de la 
mise en œuvre des matériaux à gradient fonctionnel est donné. Enfin, un résumé des travaux 
sur les méthodes de caractérisation et les méthodes d’obtention des propriétés effectives des 


matériaux à gradient fonctionnel est présenté. 


1.1. Les matériaux composites classiques 

Les matériaux composites ont été utilisés par l’homme depuis sa présence sur terre. Ils 
ne sont pas une innovation ; ils ont été utilisés dès l'antiquité par l’homme ; bois et torchis 
sont des matériaux composites de la vie habituelle. Les matériaux composites sont aussi 
utilisés pour optimiser les exploits de certaines armes ; tel que : les arcs tatars, pour lesquels 
les parties comprimées sont en corne et les parties tendus en bois collés où en tendons de 
bœuf. Les épées de Damas et les poignards Japonais dont la lame est en acier où en fer doux : 
la portion en acier est empilée comme une pâte feuilletée, menant les impuretés et les 
rétrécissements selon le sens longitudinal puis façonnée en U dans lequel est placé le fer 
doux. Dans ce cas, le sabre résiste alors à la flexion et aux chocs. Nous voyons 
subséquemment prendre forme, à cette ère, la subdivision des produits composites à grande 
diffusion et des pièces composites à hautes performances. Au sens large, le mot "Matériau 
composite" signifie ” un matériau constitué de deux ou plusieurs parts différentes”. Il est le 
résultat de l'assemblage de deux matériaux de natures distinctes, se complétant afin d'aboutir 
à un matériau dont l'ensemble des performances est supérieur à celui des composants pris 
séparément. Dans le cadre du présent travail, la nomination “matériau composite" ou 
"composite" est : un matériau composite constitué d’une phase discontinue, de renfort, rigide 
et résistante par rapport à la phase continue de la matrice. De ce fait, ses propriétés 
mécaniques découlent essentiellement de la nature des constituants et leurs propriétés, de la 
géométrie du renfort (forme, élancement) et de sa distribution (répartition, orientation), en 


conséquence de sa concentration (taux volumique ou massique). 
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1.1.1. Classement des matériaux composites 


Nous pouvons classer les matériaux composites selon de nombreuses manières. Selon 
la nature de la matrice, les matériaux composites se classent en trois catégories : 

Composites à matrice polymérique : ils ont été développé particulièrement pour les 
applications aéronautiques où la réduction de la masse est fondamentale. Les polymères sont 
donc caractérisés par une faible densité, une résistance mécanique corrélativement faible et 


une importante déformation à la rupture. 


Composites à matrice métallique : La matrice pour ces cas est métallique 
(aluminium ou titane), elle est renforcée par des renforts couramment non-métalliques, 
fréquemment des céramiques. Cette gamme de composite a des propriétés mécaniques 
meilleures ou plus adaptables au chargement que ceux à matrices polymériques. Leurs 


domaine d’application est la fabrication des moteurs d'automobile. 


Composites à matrice céramique : Pour ce type de composite la matrice est de type 
céramique telle que le verre et le carbure de silicium. La matrice est associée à des renforts 
tels que les métaux et le carbone. Le but essentiel du développement de ces matériaux 
composites est d'améliorer les propriétés mécaniques telles que la ténacité et la résistance 


au choc thermique. 


Une autre manière de catégorisation des matériaux composites se base sur le type de 
renfort figure 1.1. Dans ce cas naissent les quatre catégories suivantes : 
Composites fibreux : Les fibres de différents types et de différentes origines 
constituent des renforts pour ce type de composites. La longueur des fibres varie selon 
l’application visée en fibres continus ayant la longueur de la pièce et en fibre courte de petite 


longueur. 


Composites stratifiés : Ils sont constitués d’au moins de deux couches pouvant être 
composées de différents matériaux similaires semblablement dans les métaux plaqués ou de 
même matériau composite empilé suivant différentes orientations comme pour des 


ES 


composites stratifiés à renforts de fibres longues. Ce dernier est une classe mixte de 


composite, comportant à la fois le composite à renfort fibreux et la technique de 


stratification. 
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Composites à renforts de particules : Le renfort est dit de particule si toutes les 
dimensions du renfort sont petites devant les autres dimensions du matériau et 
approximativement égales. Les particules dures sont éparpillées aléatoirement dans la 


matrice qui est moins rigide. 


Composites à renforts de paillettes : Les dimensions des paillettes sont très faible 
par rapport aux autres dimensions. L’éparpillement de ces "particules minces” est 
généralement aléatoire. Cependant, les paillettes peuvent être ajustées parallèlement dans le 


but d'avoir des propriétés dans le plan plus homogènes. 


composite à renforts composite à renforts composite à renforts composite stratifié 


de particules de paillettes de fibres 


Figure 1.1 Classification des matériaux composites en fonction du type de renfort 


1.1.2. Constitution des matériaux composites 

Généralement un matériau composite est constitué de matrice et de renfort. Si le 
matériau composite est non endommagé, les renforts et la matrice sont absolument liés sans 
y avoir ni glissement ni séparation entre eux. Les renforts sont sous forme de fibres continues 
ou discontinues, particules où paillettes. Le rôle du renfort est d’assurer la fonction de 
résistance mécanique pour prendre les efforts extérieurs. La matrice assure la fonction de 
cohésion entre les renforts de manière à répartir les sollicitations mécaniques, tandis que 
l’arrangement des fibres et leur orientation permettent de renforcer les propriétés mécaniques 
de la structure. Les propriétés mécaniques à l’interface entre fibres et matrice sont de grandes 
importances dans l’exécution d’un composite. Pour obtenir de bonnes caractéristiques 
mécaniques, 1l ne doit y avoir n1 glissement n1 séparation entre les différents constituants du 


composite. 


1.1.2.1. Les matrices 


Dans la majorité des cas, la matrice composant le matériau composite est une résine 
polymère. Les résines polymériques existent en grand nombre et chacune à un domaine 
particulier d'application. Lorsqu’une tenue de la structure aux très hautes températures est 


indispensable, des matériaux composites à matrice métallique, céramique ou carbone sont 


15 


Chapitre 1 Etude bibliographique sur les matériaux composites 


préférables. Une classification, des types de matrices, couramment rencontrées dans la 


littérature se présente comme suit : 


Résines thermodurcissables 
Ces résines sont les plus utilisées de toutes les matrices, pour des applications 
commerciales, industrielles et de transport. Elles ont une excellente adhésion avec les 


différents fibres et particulièrement ceux en verre. 


Résines thermoplastiques 

S1 la résine thermodurcissable est soumise à des températures élevées, il se génère des 
relations chimiques au cours desquelles des liaisons covalentes sont formées entre les 
chaînes des molécules résultantes. Le polymère final se présente sous la forme d'un corps 
solide infusible. La polymérisation est irréversible, c’est à dire que le matériau ne peut être 


mis en forme qu'une seule fois, au moment de sa fabrication. 


1.1.2.2. Les fibres 


La rupture dans les matériaux à haute résistance ou à haut module d’élasticité est 
fréquemment provoquée par la propagation de fissures. Les matériaux fibreux sont 
fondamentalement plus résistants à la rupture que ceux en forme massive, puisque la taille 
des défauts est initialement limitée par le diamètre faible. Pour un composite fibreux, la tenue 
mécanique est assurée particulièrement par les fibres. Suite à sa nature filamenteuse, la 
rupture de quelques fibres a pour résultat le transfert du chargement sur les autres fibres, ce 
qui prévient la rupture catastrophique de la structure. Les fibres les plus souvent rencontrées 


dans les composites sont les suivantes : 


Les fibres de verre : Les fibres de verre sont très répandues dans les applications à 
basse où à haute performance telles que les réservoirs de propulseurs de fusée. Leurs 
avantages impliquent un prix compétitif, une disponibilité et une résistance élevée. 
Néanmoins, à cause de leur rigidité corrélativement faible, les fibres de verre sont 
progressivement remplacées par des fibres d’aramide ou des fibres de carbone dans les 
applications de haute performance. 

Les fibres de carbone : Maintenant, les fibres de carbone sont les plus répandues 
pour les matériaux composites à hautes performances. Deux avantages principaux jouent en 
faveur de ces fibres : leur fabrication plus adaptée à la production à grande échelle que 


d'autres fibres à hautes performances et leurs excellentes propriétés mécaniques plus 
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facilement transférables aux matériaux composites. Leur prix reste toutefois excessif pour 
les produits à grand public. 

Les fibres d’aramides : Ces fibres appartiennent à la catégorie des fibres polyamides 
aromatiques. Elles possèdent une résistance élevée et une rigidité amplement supérieure à 
celle des fibres de verre. Leur latitude aux dommages est également très bonne. Leur 
désavantage principal est leur résistance en compression inférieure à celle des fibres de 
carbone et une cohésion relativement faible aux matrices. 

1.1.3 Réflexions sur l’usage des matériaux composites 

Dans la conception des produits finaux en matériaux composites, 1l est indispensable 
d'évaluer et de comparer les composites en matériaux conventionnels pour bien adopter les 
matériaux à utiliser. Les avantages et les inconvénients principaux des matériaux composites 
sont présentés ci-dessous. 

1.1.3.1. Avantages des matériaux composites 


Les avantages les plus cités des matériaux composites sont : 


e Adaptation de leurs propriétés mécaniques 
Afin d'obtenir un composite capable de résister à diverses sollicitations, 11 faut alors 
associer différentes orientations de fibres. Un stratifié est un empilement de plis élémentaires 
se comportant ensemble comme un élément structural. Un pli élémentaire est anisotrope, ou 
orthotrope dans la majorité des cas, avec une résistance et une rigidité dans la direction des 
fibres beaucoup plus élevées que dans l’autre direction. L’avantage capital du stratifié 
composite est que les plis élémentaires peuvent être orientés de telle façon que la résistance 


dans une direction donnée concorde avec les chargements prévus. 


e Haut module d’élasticité et haute résistance mécanique 
Les propriétés mécaniques élevées, surtout la résistance et la rigidité des matériaux 
composites à hautes performances, permettent de répondre aux exigences de l'industrie. 
D'autres industries en profitent également, par exemple la vitesse d'une balle de tennis est 


plus élevée avec des raquettes en carbone/époxy. 


+ La légèreté 
Les matériaux composites à haute performance ont une densité de l'ordre de 1,6, en 
revanche celle des alliages d'aluminium est de 2,7. Ils autorisent donc un gain de 10 à 30% 
sur les éléments de structure. Ce bénéfice de masse permet d'employer des moteurs moins 


puissants, qui consomment moins et permettent de réduire la taille des réservoirs de 
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combustible pour le même cahier des charges de l'avion. Il permet aussi des rayons d’action 
plus grands. 
e Facilité de Production 
Les matériaux composites permettent la simplification de l'assemblage dans les 
structures, ceci corrige relativement leur prix élevé. La diminution du nombre de pièces par 
rapport aux matériaux ordinaires peut être consistante. 
e Facilité de Maintenance 
Les composites nécessitent moins d'entretien que les alliages métalliques. D'une part, 


ils ne sont pas sensibles à la corrosion et d'autre part leur tenue en fatigue est très bonne. 


1.1.3.2. Les inconvénients des matériaux composites 


Les matériaux composites ne sont pas une solution miracle pour toutes les applications 
malgré leurs avantages impressionnants. Des inconvénients et/ou des problèmes existent et 


peuvent empêcher leur usage. Les inconvénients les plus usuels sont : 
Prix de revient élevé (Coût élevé) 


Les matériaux composites ont été développés particulièrement pour répondre aux 
besoins de l’industrie aérospatiale qui est très exigeante vis-à-vis de la sécurité et de la 
technologie utilisée, ce qui rend lent l'élargissement de l’utilisation des composites aux 
produits plus populaires (de grande consommation) à cause de leur coût élevé. 

Conception et analyse sophistiquées 

A la différence des matériaux classiques, profondément étudiés, les matériaux 
composites sont fréquemment à la fois hétérogènes et anisotropes ; rendant ainsi difficile 
leur conception et leur analyse. C'est pourquoi, ils nécessitent de nouvelles approches de 
conception, d’étude et d’analyse, généralement plus approfondies et plus sophistiquées. 
L'hétérogénéité exige l'analyse de ces matériaux de deux points de vue, la première 
macroscopique et la deuxième microscopique : la macro-mécanique et la micromécanique. 
De même, les matériaux anisotropes requièrent plus de propriétés mécaniques que les 
matériaux isotropes pour établir les relations contraintes-déformations. Dans le cas des 
composites stratifiés, la rupture se produit de plusieurs façons différentes. Il est donc difficile 
de mettre une théorie de résistance tenant compte de tous les modes de rupture et de leurs 
interactions. Par conséquent, la détermination des résistances est plus complexe nécessitant 


ainsi des vérifications expérimentales souvent onéreuses. De plus, la résistance d'un stratifié 
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peut être influencée par les contraintes de cisaillement transversal entre les couches (inter- 
laminaires) au voisinage des bords libres (problème de délaminage). 
e  Fragilité et Faiblesse des Assemblages 

L'assemblage des pièces en matériaux composites pose toujours des problèmes de 
faible résistance, à cause de la concentration des contraintes, et cela vu que les matériaux 
composites sont plus fragiles que les matériaux métalliques. Par conséquent, Les 
assemblages sont souvent renforcés par des inserts métalliques ou par l'augmentation de 
l'épaisseur du composite au droit des singularités entrainant ainsi une masse supplémentaire 
sur la structure. 

e Latitude aux dommages 

Un des points faibles les plus sérieux des matériaux composites est la latitude aux 
dommages. Des dommages inévitables de natures diverses peuvent être engendré durant la 
vie d'une structure soit en service ou pendant la maintenance. D’une manière générale, plus 
un matériau est ductile, plus il est capable de tolérer un dommage, vu que la ductilité a une 
capacité d'absorber de l'énergie. Par conséquent, sous un choc, les structures métalliques ont 
une aptitude de se déformer plutôt que de se briser. Le caractère fragile des matériaux 
composites n’admet pas, en revanche, de subir un choc sans subir de dommages. Les 
atteintes sont souvent des fissurations internes de la matrice, indétectables sur la surface de 
la structure. Ce type de dommages sur la matrice réduit amplement la résistance en 
compression des pièces endommagées. Les atteintes des fibres diminuent la résistance en 
compression et en traction du composite. Un programme d’inspection (analyse de la 
résistance résiduelle et analyse de la propagation des dommages) est nécessaire pour évaluer 


la tolérance aux dommages des structures en matériaux composites. 


1.1.4. Les mécanismes de rupture dans les composites 

Le matériau composite à renforts de fibres longues est très coutumier dans les 
applications où la réduction de masse est cruciale. Plus l'utilisation d’un matériau est élargie, 
plus la probabilité de rupture éventuelle est élevée. La capacité à déterminer le type des 
ruptures, par exemple en termes d’identification des modes de rupture, des paramètres 
caractéristiques, ou des valeurs critiques à la rupture, est indispensable pour assurer 


l'intégrité des pièces en service et pour la conception des futurs produits. 


1.1.4.1. Modes de rupture dans les composites unidirectionnels 
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La rupture dans les composites stratifiés se produit de différentes façons très 
complexes. Les modes de rupture découlent de la stratification et de la direction du 


chargement par rapport à l'orientation des fibres. 


Les figures 1.2 et 1.3 présentent les allures des ruptures par flexion et par traction 
respectivement. 


PR— œ) 
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Figure 1.3 Ruptures par traction de différentes séquences d’empilement (0°-90°) 


La description de la rupture à l'échelle du pli sert au classement des mécanismes de 
rupture. Les stratifiés à renforts de fibres longues montrent trois types de rupture : rupture 
intra-laminaire, rupture inter-laminaire et rupture translaminaire. Ces trois modes de rupture 
sont présentés à la figure 1.4. Elles définissent le plan de rupture par rapport aux composants 
du matériau. La rupture intra-laminaire se manifeste à l'intérieur d'un pli tandis que la rupture 
inter-laminaire est une rupture entre deux plis contigus. La rupture translaminaire a une 


orientation transversalement à l'orientation des fibres pour le pli endommagé. 
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Rupture intralaminaire Rupture interlaminaire Rupture translaminaire 


Figure 1.4 Les trois modes de rupture dans un stratifié 0°/90°/0° (Smith 1993). 


e La rupture intra-laminaire : Elle est due particulièrement à la faible résistance 
de la matrice et/ou de l'adhérence entre la matrice et les fibres. Elle est provoquée par les 
contraintes se trouvant dans le plan du stratifié. Un pli se dégrade par la résultante de la 
contrainte de traction dans la direction normale aux fibres. Ce type de rupture est connu sous 
l’appellation de "fissure transverse". Dans la majorité des cas, ce type de fissure de la 
matrice se manifeste avant la rupture des fibres. 

e La rupture inter-laminaire se manifeste à l'interface entre deux plis d'un 
empilement. Ce type de rupture se caractérise par la rupture de la matrice et la décohésion 
fibre/matrice. Ce mécanisme de rupture inter laminaire contient peu de rupture de fibres. De 
même que pour le cas des matériaux métalliques, la rupture inter laminaire se produit selon 
des modes distincts connus en mécanique de la rupture ; à savoir le mode I (ouverture), le 
mode II (cisaillement), le mode III (cisaillement anti-plan) ou une combinaison des trois 


modes. 


MODE 1 MODE I MODE IN 
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Figure 1.5 Les trois modes principaux de rupture dans les composites 


e La rupture translaminaire : concerne la rupture de fibres. Les surfaces de 
rupture sont donc généralement marquées par la morphologie rugueuse des bouts de fibres. 
En effet, la contrainte de rupture des fibres est plus importante que celle de tous les autres 
constituants d'un composite stratifié. Par conséquent, ce mécanisme de rupture entraîne 


souvent la rupture totale du stratifié. La rupture translaminaire peut être séparée en deux 
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modes selon les chargements : la rupture en traction et le micro-flambage en compression. 
La rupture peut être provoquée par un mode individuel ou une combinaison de deux modes. 

1.1.4.2. Risque du phénomène de délaminage 

L'un des avantages essentiels des composites stratifiés à renforts de fibres longues est 
l’aptitude d'orienter les fibres de chacun des plis afin d’obtenir des propriétés, souvent la 
résistance et la rigidité, adaptées aux chargements dans les directions attendues. Malgré les 
excellentes propriétés dans le plan, les composites stratifiés présentent un problème propre 
aux matériaux réalisés par empilement (stratification) ; il s’agît de la rupture inter-laminaire. 
Ce mécanisme de rupture connu sous l’appellation de "délaminage" se caractérise par un 
décollement et/ou une décohésion entre les couches du stratifié. Un stratifié sollicité par un 
chargement présente différents taux de dégradations. Dans le schéma d'évolution des défauts 
les plus “classiques”, la matrice et l'interface fibre/matrice sont les premières à s’abîmer 
(Figure 1.6-a). Les premières anomalies sont en conséquence la microfissuration de la 
matrice et la décohésion fibre/matrice à l'échelle microscopique. Puis, ces défauts 
s'accroissent de manière permanente à l'échelle du pli par coalescence (Figure 1.6-b), les 
micros dégradations s’accumulent pour constituer des fissurations transverses apparentes. 
Ces derniers (les fissures transverses) peuvent arriver à l'interface du pli et provoque le 
délaminage sous l'effet des contraintes inter-laminaires (Figure 1.6-c). Généralement, le 
délaminage et son évolution dépendent du nombre de plis, de leur empilement et du 


chargement. 


“(e) Micro-défauts | : (b) Coalescence (ce) Fissuration transverse 
Figure 1.6 Développement et Evolution des défauts 


Fissuration 


Couche à 90° 


Délaminage 


Fissuration 


[a ouche à Q° \ 
Rupture des fibres longitudinale 


Figure 1.7 Le délaminage dans un composite 
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Les principaux micro-mécanismes d’endommagement accompagnant un délaminage 
se présentent schématiquement comme dans la Figure 1.8. Ils contiennent : 

e Une zone endommagée : La forte concentration de contraintes tout autour de 
la pointe de fissure fait naitre une zone de dommage où se trouve la déformation plastique 
et/ou des microfissures dans la matrice. 

e Une fissure latérale : Les microfissures dans la zone endommagée peuvent se 
transformer en des fissures frontales dans la matrice aux alentours du plan de délaminage. 

e Le pontage de fibres : L’existence de fissures en dessous et/ou au-dessus du 
plan de délaminage facilite la création de passerelles de fibres reliant les deux surfaces 
délaminées. Notons que certains ponts de fibres se rompent pendant l'avancée du 
délaminage. 


Side 
Cracks 
+ Length }— 
£ 


|— 
FE 


Profile | \ 
Polymer Crack Bridging  Plastically Deformed/ 
Matrix of Fibers Microcracked Matrix 


Figure 1.8 Les mécanismes d’endommagement accompagnant le délaminage 


1.2. Les matériaux à gradient fonctionnel «FGM> 

Pendant un certain temps l’industrie aéronautique, navale et autres applications 
d’ingénierie ont utilisé des matériaux composites légers avec un rapport élevé de résistance- 
poids. La figure 1.9 représente les différentes catégories de matériaux en fonction de leurs 
densités p et de leurs modules d’élasticité (module de Young) E. Cependant, le matériau 
composite traditionnel (classique) est inapte d'être employé dans des environnements à haute 
température. Tandis que les matériaux métalliques ont été utilisés dans le domaine industriel 
durant de nombreuses années en raison de leur bonne résistance et de leur ténacité. De même 
que pour un matériau composite, à température élevé, la résistance du métal se réduit 
pareillement à un matériau composite traditionnel (classique). Contrairement aux métaux et 
aux matériaux composites classique, les matériaux céramiques ont d'excellentes 
caractéristiques de résistance à la chaleur (haute température). Cependant, l’utilisation des 


céramiques est souvent limitée en raison de leur faible résistance : mécanique et au choc. 


Récemment, une classe originale de matériaux composites nommée : Matériaux à 


Gradient Fonctionnel (MGF) (FGM en langue anglaise) a suscité une réflexion considérable 
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des chercheurs dans le domaine des composites. Un matériau composite de type MGF est un 
composite non homogène fabriqué à partir de différents composants de matériaux 
(généralement du céramique et métal), avec un effet considérable de couplage membrane- 
flexion du point de vue comportement mécanique. La figure 1.10 illustre un exemple de ce 
type de matériau tiré de Yin 2004. Elle montre l’intégration de particules de différentes 


formes (généralement sphérique) dans une matrice isotrope. 


Céramiques 
poreux Métaux et 
alliages 


Polymères 


p (kg/m°) 


Figure 1.9 Classification des matériaux en fonction de leurs densités p et de leurs 
modules de Young E (Ashby 2005) 
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Figure 1.10 Matériau à gradient fonctionnel (Yin 2004) 


Une micrographie par microscope à balayage électronique d’une section transversale 


d’un matériau FGM est présentée sur la figure 1.11. 


" 
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ALO3 < > SUS304 


Figure 1.11 Micrographie par microscope électronique à balayage d'une section 
transversale d’un fonctionnel (Al1O3-SUS304). 


Pour un matériau à gradient fonctionnel (FGM) les différentes phases micro- 
structurelles ont des fonctions différentes. Le matériau à gradient fonctionnel absolu atteint 
l’état multi structural par échelonnement de leurs propriétés. En variant graduellement la 
fraction volumique des constituants du matériau, leurs propriétés matérielles montrent un 
passage lisse et continu d’une surface à l’autre, éliminant ainsi les problèmes d'interface 
(inter-laminaire) et l'atténuation des concentrations de contraintes thermiques. Cela du fait 
que le constituant céramique du matériau à gradient fonctionnel est capable de résister à des 
hautes températures, tandis que le constituant métallique garantis une meilleure performance 
mécanique et réduit la potentialité de la rupture brutale catastrophique. Donc, en utilisant les 
FGM, on anticipe que les propriétés du composé vont être améliorées et de nouvelles 
fonctionnalités vont voir le jour. Les premiers matériaux à gradient fonctionnel ont été 
réalisés par des chercheurs Japonais au milieu des années quatre-vingt (Yamanoushi 1990 et 
Koïizumi 1993). Depuis, des efforts considérables ont été fourni pour développer des 
matériaux à haute résistance basés sur les matériaux à gradient fonctionnel. Ces produits ont 
servie initialement comme des barrières thermiques pour les structures aérospatiales et les 
réacteurs de fusion ; Hirai and Chen 1999, Chan 2001 et Uemura 2003. Actuellement, on 
les utilise pour un usage général en tant que composants structuraux dans un environnement 
de haute température. Les applications éventuelles des matériaux à gradient fonctionnel sont 
variées et nombreuses Müller 2003, Qiu et al. 2003, Liu 2003, Watari et al. 2004 et 
Paszkiewicz et al. 2008. 


Plusieurs articles dans des revues de recherche traitant divers aspects sur les FGM ont 
été publiés au cours des dernières années (Fuchiyama 1995, Markworth 1995, Tanigawa 
1995, Noda 1999 et Paulino 2003). Ils affirment que la plupart des recherches sur les FGM 


avaient davantage porté sur l'analyse des contraintes thermiques et mécaniques de la rupture. 
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Une étude complète menée en flexion, flambage et l'analyse des vibrations des plaques et 
coques en FGM a été présentée par Shen 2004. En 2007, Birman et Byrd ont fait une synthèse 
de recherche bibliographique sur les principales évolutions en FGM incluant le transfert de 
chaleur, contrainte, stabilité et analyses dynamiques, expérimentation, fabrication et 


conception ainsi que les applications et rupture. 


1.2.1 Notions sur le gradient des FGM 

À l'inverse aux matériaux homogènes, les propriétés des FGM changent non seulement 
en fonction de leur composition mais également en fonction de la connectivité du réseau 
interne de la structure. Sur la figure 1.12, nous montrons une représentation schématique de 
la microstructure d’un matériau FGM suite à l’augmentation du contenu de la deuxième 
phase. Pour des fractions faibles, la deuxième phase existe sous forme de grains isolés 
éparpillées dans une matrice (a). Quand le contenu de la deuxième phase accroit, les 


particules débutent à être en contact et à former des tas agglomérés (b). 


Si l’augmentation de la phase deux continue, une transition microstructurale critique 
aura lieu, dans ce cas la deuxième phase n’est plus éparpillée, mais elle devient plutôt reliée 
entre elles sur des volumes importants (c, d, e). Cette phase de transition a un effet important 
sur l’évolution des propriétés du matériau FGM final. Nous remarquons qu’une petite 
transformation de la composition du FGM va se traduire par une variation distincte des 


propriétés telle que la conductivité thermique ou électrique (Morinaga 1993). 


Figure 1.12 Changement de la microstructure d’un FGM 


Les FGM se caractérisent par les trois types de gradients suivants : 
Le Gradient de la porosité (Taille et pourcentage des pores) 


Le gradient de la porosité s’obtient soit par le dépôt du mélange des poudres à 
différentes formes des particules ou par variation des paramètres de dépôt en utilisant des 
détenteurs d'espace. La taille des pores (vides) se crée en général en variant la taille des 
particules. Les différents comportements au frittage des poudres peuvent être utilisés pour 


produire un gradient de porosité. 
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Le Gradient de composition chimique d’un matériau monophasé 
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Le dépôt des poudres avec un changement continu de la composition de la mixtion 
durant le frittage entraine la naissance d’une matière à phase unique avec un changement 
flexible de la répartition. 


Les gradients de la teneur en volume des phases et les gradients de taille de grain 
(Matériaux multiphasique) 

La majorité de la bibliographie se centralise sur le gradient en volume des contenants 
des distinctes phases. La microstructure résultante du FGM est constituée de deux ou 
plusieurs phases avec des gradients en volume des constituants. La combinaison des phases 
peut admettre des systèmes : métal-métal, métal-céramique et céramique-céramique. 


1.2.2 Les méthodes de caractérisation des FGM 


À l'égard de la généralisation des applications progressivement croissantes des FGM, 
il est nécessaire de trouver une méthode adéquate et fiable pour quantifier et mesurer les 
propriétés physique et mécanique des matériaux FGM en utilisant des essais non destructifs. 
Il est donc important de développer et d'avoir une technique fiable et efficace de 
détermination des propriétés des différents matériaux entrants dans la conception des FGM 
et 1l est particulièrement important d'avoir une technique pour évaluer les FGM après leur 
fabrication et leur mise en service. Ceci nous permet de vérifier si la propriété réelle du 
matériau concorde à l'exigence de conception pour satisfaire la sécurité et la qualité. En 
revanche la microstructure d'un produit FGM varie dans l’espace, par conséquent, la ténacité 
et la résistance d'un FGM ne peuvent être quantifiées naturellement par les tests classiques 
de la mécanique et la mécanique de la rupture ordinaire. Afin de déterminer le module 
d’élasticité (module de Young) de l'assortiment d'un FGM, habituellement à partir de la 
courbe contrainte-déformation, une jauge de déformation peut être fixée à la surface de 
l'échantillon au cours d'un essai de flexion à quatre points. La corrélation du module de 
Young de la combinaison d'un FGM s’estime en mesurant les fréquences de résonance de 
flexion d'une éprouvette rectangulaire à l’aide de la technique de la résonance forcée. Cette 
dernière, la résonance forcée, s’articule sur l’utilisation des champs d'ondes élastiques dans 
les structures et les procédures de calcul inverse qui sont encourageante pour la 
caractérisation des FGM ; puisque aujourd'hui des solveurs puissants sont à la portée des 
chercheurs pour l'analyse des ondes élastiques se propageant dans des milieux élastique telle 
que les FGM. 


1.2.3. Design, optimisation et conception inverse des FGM 


La conception d'un FGM est basée généralement sur les conditions de leur service. De 
ce fait, des efforts importants sont nécessaires pour sélectionner une combinaison de 
matériaux optimale dans le but d’alléger une contrainte (thermique par exemple) où une 
distribution optimale inter-matériau de la composition convenant aux conditions de 
production. Généralement, le concepteur est confronté à concevoir un système pour une 
performance optimale donnée sous certaines conditions d’usage et de contraintes 
d’utilisation et davantage, il a un éventail de choix du matériau à utiliser. Dans le cas des 
composites stratifiés fibreux, il a l’opportunité additionnelle d'orientation des fibres d’un 
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stratifié en fonction des exigences de conception. À l'égard de l'apparition des FGM, il peut 
aller encore loin en utilisant la notion de la conception où analyse inverse. 


1.2.4. Les propriétés graduelles des FGM 


Les variations graduelles des propriétés, d’un matériau FGM dans l’espace, comme 
elles ont été évoquées jusqu'à présent, sont initialement représentées par des polynômes, des 
fonctions exponentielles etc. En 1991, Fakui a étudié les effets de la fraction volumique de 
la poudre et de la force centrifuge sur le profil ayant une fonction de distribution polynomiale 
de degré trois (cubique). Puis, en 1993, Jin et Noda utilisent une variation exponentielle du 
profil de la distribution et déterminent quelques paramètres en liaison avec la fonction 
exponentielle employée. Dans la même année (1993), Matsuzaki et al. ont optimisé 
numériquement l'épaisseur et la distribution de composition d’un FGM vis-à-vis du transfert 
thermique dans des systèmes de réfrigération. De même, en 1993, Williamson d’une part et 
Drake d’autre part ont analysé les caractéristiques d’une microstructure en considérant les 
contraintes résiduelles. 


1.2.5 Relations constitutives et modélisation du gradient 


Généralement, deux façons de modéliser la graduation des propriétés de matériau dans 
les solides sont à considérer : 


(1) L’approche macroscopique, elle suppose un profil de variation des fractions 
volumiques 
(11) L’approche micromécanique pour l’étude des milieux non homogène. 


Pour la première approche (macroscopique), des fonctions polynomiales de différents 
degrés sont utilisées. Mais à l’échelle microstructural, le matériau à gradient fonctionnel 
se distingue par une transition à partir d'une phase dispersée à une autre phase à réseau 
inclusivement. Il est à noter que durant la dernière décennie, la micromécanique des 
milieux aléatoires inhomogènes a fait l'objet de plusieurs études intensives. Toutefois, des 
études semblables, ont principalement tenu sur des inclusions ou l'inhomogénéité de 
différentes formes du substrat homogène (homogénéité statistique) est considérée. 


Propriétés matérielles 


Figure 1.13 Illustration des propriétés matérielles des profils. 
1.2.6 Domaines d’applications des FGM 


Les FGM se caractérisent par un large domaine d’application. Sur le tableau 1.1 nous 
présentons les variétés des domaines d'applications réelles et éventuelles des FGM dans : les 
outils de coupe, les systèmes de transport, l’optique, pièces de machines, semi-conducteurs, 
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et le biosystème. Les applications éventuelles abritent tous les cas pour lesquels une 
combinaison de fonctions incompatibles est indispensable. 


Tableau 1.1 Les divers domaines d’applications des FGM. 


Conversion d'énergie | Générateur thermoélectrique, pile à combustible, une cellule 


solaire 
Matériaux de Blocs en béton résistant au feu, rayonnement, absorbant de 
construction l’humidité, milieux agressive. 
Ingénierie Outil coupant, les arbres, rouleau, lame de turbine, Matériaux de 
Construction 
Énergie Nucléaire Composants de réacteurs nucléaires, première paroi de du 


réacteur fusion 


Aéronautiques Composants de véhicules spatiaux, corps d’avion spatial 

Aérospatiale 

Biomatériaux Implant, peau artificielle, le système de délivrance de 
médicaments 


Usines de produits Ç | 
EÉchangeur de chaleur, cuve de réaction 


chimiques 

optique Fibres optiques, lentilles 

Électronique Semi-conducteur de bande graduée, capteur 
Produits de base Articles de sport, carrosseries, vitres 


1.2.7 Les propriétés effectives des matériaux à gradient fonctionnel 

La plus part des FGM sont fabriqués à partir de deux phases (types) de matériaux ayant 
des propriétés distinctes. Un aperçu détaillé d’une microstructure échelonnée tangible est 
ordinairement non disponible, sauf éventuellement pour des renseignements sur la partition 
et la distribution de la fraction volumique. Comme la fraction volumique de chaque phase 
change et évolue graduellement suivant l’orientation de la progression, les propriétés 
effectives des FGM changent tout le long de cette direction. À cet effet, il existe deux 
approches pour construire des modèles types de matériaux à gradient fonctionnel. La 
première, consiste à prendre une variation par morceaux (variation en escalier) de la fraction 
volumique que ce soit de la céramique et/ou du métal. Dans ce cas le FGM garde la même 
fraction volumique pour chacune des régions (des couches quasi-cohérentes de métal- 


céramique Figure 1.14). La deuxième démarche, sert à prendre une variation continue de la 
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fraction volumique des deux composants, l’une à l’inverse de l’autre, semblablement à une 
fonction dépendant de la coordonnée suivant laquelle varie la propriété considérée 


(généralement suivant l’épaisseur z). 


“(y = 


Avec h l’épaisseur de la structure (poutre, plaque, coque) et p (0 < p < «) est un paramètre 


lié à l’exposant de la fraction volumique et qui gère le profil de la variabilité du matériau le 
long de l’épaisseur de la strate du FGM. Il existe un nombre infini de la distribution de la 


composition en faisant varier le paramètre de l’exposant p (figure 1.15). 


(a) (b) 


Figure 1.14 Modèle analytique de distribution dans une couche de matériau FGM 


Figure 1.15 Variation de la fraction volumique du métal à travers l'épaisseur. 
Dans le but de faire une analyse approfondie et une modélisation exacte des structures en 


FGM leurs propriétés doivent être représentées avec justesse. Cette dernière s’obtient par la 
simple règle de mélange, très connue dans le domaine des matériaux composites, connu sous 


le nom de ‘modèle de Voigt’. Les propriétés des matériaux effectives P- d’un FGM, comme 
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le module d’élasticité E;(module de Young) et le coefficient de dilatation thermique 


af S’expriment par la formule (1.2), suivante : 


Pr — P;Vr; (1.2) 
i=1 


Avec P; et Vr; sont les propriétés du matériau et la fraction volumique du matériau 


constitutif i respectivement. 


La sommation des fractions volumiques de tous les constituants du matériau FGM est 


> Pei (1.3) 
i=1 


Les structures en FGM sont conçues pour être utilisé dans des environnements à haute 


égale à l’unité. 


température, pour lesquels les propriétés mécaniques des matériaux de base constitutifs des 
FGM ne sont plus permanentes, mais elles prennent des variations importantes en fonction 


de la température (Reddy 1998). 


Il est donc important de considérer amplement cette relation de dépendance de la 
température pour une prédiction précise de la réponse mécanique. De ce fait, le module de 
Young Ef, le coefficient de Poisson vf, le coefficient de dilatation thermique as, et la 
conductivité thermique Kkf sont donc fonction de la température. 


Tableau 1.2 Propriétés à température ambiantes 


Matériaux Propriétés mécaniques 
E [GPa] V p[Kg/m°] 
Aluminium (Al) 70 0.3 2702 
Alumina(Al, 03) 380 0.3 3800 
Zircon (ZrO;) 200 0.3 5700 


Dans le cas d’une distribution en loi de puissance suivant la direction de l’épaisseur de 
la structure (poutre, plaque où coque) les propriétés matérielles effectives E;, v;, &; et K; 
s’évaluent par les expressions (1.4) jusqu’à (1.7) respectivement (Reddy 1998). 


2z +h\? 1.4 
B= TE ECS) + EC) De 


2z + h\? 15 
= But = IE) + 00) . 


2z+h 


p 1.6 
pr = Len PI) + re) dE 
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ap = Lan() - &(TI( SE) + ae(n) ss 


Il est clair que les propriétés Er, pr, et vf sont fonction de la position et de la 


température à la fois. 


Cette représentation est simple et facile à appliquer pour prédire les propriétés des 
matériaux et les réponses globales, mais, en raison de ces simplifications la conformité est 
affectée par la microstructure graduée. Comme il est mentionné précédemment les 
renseignements précises sur la forme, la taille et la répartition des particules ne sont pas 
disponibles et les modules élastiques effectifs des microstructures graduées doivent être 
évalués en fonction de la distribution de la fraction volumique et la forme approximative de 
la phase éparpillée. Des modèles micromécaniques ont été développés cette dernière 


décennie dans le but d’évaluer les propriétés effectives des FGM. 


La théorie de Mori-Tanaka (Mori et Tanaka , 1973)pour évaluer les propriétés effectifs 
est approprié pour les zones de la microstructure à gradient, ayant une matrice continue 
convenablement définie et une phase discontinue telle que représentée à la Figure 1.10. Cette 
théorie de Mori-Tanaka prend en compte l'interaction des champs élastiques entre inclusions 
voisines. L’approche auto-consistante (self-consistent method) (Hill, 1956) se base sur le 
faite que chaque inclusion de renforcement est submergée dans un milieu continu dont les 
propriétés effectives sont celles du composite. L’approche en question ne fait pas de 
distinction entre la matrice et les inclusions (phases) de renforcement. Les modules globaux 
sont espérés dans un autre composite dans lequel les rôles des phases sont commutables. 
Cela rend le matériau particulièrement approprié pour déterminer les modules effectifs dans 


les régions ayant une microstructure interconnectée du squelette Figure 1.16. 


Figure 1.16 Squelette Microstructure du matériau FGM (Vel et Batra 2002) 


Ainsi, le modèle de Mori-Tanaka témoigne des prévisions précises des propriétés avec une 


matrice continue convenablement définie et des inclusions discontinues. Tandis que le 
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modèle auto-consistant est mieux adapter pour les microstructures squelette caractérisées par 


une large zone de transition. 


1.3 Conclusion 

Dans le présent chapitre nous avons fait un tour d’horizon sur les matériaux composites 
et plus précisément sur les matériaux composites à gradient fonctionnel. Dans le chapitre 
subséquent nous abordons la modélisation des structures poutres en FGM par la méthode 


des éléments finis. 
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ANNEXE 


Propriétés matérielles et lois de distribution des structures en FGM 


1. Propriété matérielle d’une structure P-FGM (Loi de Puissance) 
La fraction volumique de la distribution obéissant à une distribution en loi de puissance 


P-FGM s'écrit comme suit : 


z+h/2\? (1.8) 
(#2) 
Où p est un paramètre matériel et h est l’épaisseur de l’élément de structure. Une fois la 
fraction volumique locale V.(z) est déterminée, les propriétés matérielles de la loi P-FGM 
se déterminent par la loi des mélanges comme suit : 
E(z) = (Ec — En)VC) + En 
z+h/2\? (1.9) 
E(z) = (Ec — En) TE En 
Où E, et E} sont les modules de Young de la surface supérieure (z = +h/2) et de la surface 


inférieure (z = —h/2) respectivement. 


n - 
00 02 04 06 08 LO 
Fraction volumique Ve 


Figure 1.17 Variation de la fraction volumique à travers l’épaisseur selon la loi P-FGM 


La variation de la fraction volumique à travers l’épaisseur d’un élément de structure P- 
FGM est représentée sur la (Figure 1.17). Il est claire que cette variation évolue rapidement 
près de la surface inférieure pour p < 1 et accroît rapidement près de la surface supérieure 


pour p > 1. 


2. Propriété matérielle d’une structure S-FGM (Loi Sigmoïde) 
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Pour assurer une bonne distribution des contraintes entre toutes les interfaces, Chi et 
al. en 2003 ont utilisé deux fonctions de loi de puissance, pour définir la fraction de volume, 
chacune sur un demi plan de l’élément de structure selon la profondeur z. Ces deux fonctions 


sont définies comme suit : 


p 
(2) — (=) pour :—h/2 < Z < 0 (1.10a) 
_z\P 
VC) = 1 -1(ue) nr (1.10b) 


Par la loi de Voigt le module de Young de l’élément S-FGM s’obtient par : 


E(z) = V.G)E, + [1 -V(2)]E, Pour :—-h/2<z<0 (1.11a) 
E(z) = V,(G)E, + [1-V,(2)lEn Pour :0 <z<h/2 (1.11b) 


L’allure de ces deux fonctions de puissance est sous forme de sigmoïde, de ce fait la 
loi de distribution est notée par S-FGM. Elle est représentée sur la figure 1.18. On remarque 


que cette distribution est antisymétrique par rapport au plan moyen. 


zh 


Fraction volumique V 


Figure 1.18 Variation de la fraction volumique à travers l’épaisseur selon la loi S-FGM 
3. Propriété matérielle d’une structure E-FGM (Loi Exponentielle) 


Une fonction exponentielle a été utilisée par Delale en 1983 pour représenter les propriétés 


matérielles des FGM. Elle est de la forme : 


E(z) = E,ePGth/2) (1.12a) 
pe () (1.12b) 
he 


Pour le cas du module de Young elle se représente comme suit : 
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E(e) 


05 04 “03 02 01 00 O1 O2 03 04 05 
zh 


Figure 1.19 Variation du module de Young à travers l’épaisseur selon la loi E-FGM 
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Dans le présent chapitre, nous proposons une théorie de poutre à déformation de 
cisaillement d'ordre supérieur pour la flexion de poutre fonctionnellement graduée. Nous 
utilisons une nouvelle fonction de cisaillement polynomiale qui satisfait aux conditions aux 
limites de contraintes de cisaillement nulles (conditions aux limites exactes sur la contrainte de 
cisaillement) à la fois sur les surfaces supérieure et inférieure de la poutre. Par conséquent, le 
facteur de correction de cisaillement n'est pas nécessaire. De plus, la présente théorie montre 
de fortes similitudes avec la théorie des poutres de Timoshenko dans certains concepts tels que 
les équations du mouvement, les conditions aux limites et les expressions résultantes des 
contraintes. Les équations gouvernant le mouvement et les conditions aux limites se dérivent 
du principe du minimum de l'énergie potentielle totale. Les poutres en matériau à gradient 
fonctionnel (FGM) présentent une variation régulière des propriétés du matériau en raison du 
changement continu (ininterrompu) des détails de la microstructure. La variation des propriétés 
du matériau se situe le long de l'épaisseur de la poutre et est supposée suivre une loi de puissance 
de la fraction volumique des constituants. Les solutions numériques obtenues par éléments finis 
avec la nouvelle fonction polynomiale de déformation de cisaillement sont présentées et les 
résultats obtenus ont été évalués et comparés avec les solutions existantes pour vérifier la 
validité de la présente formulation. Enfin, les influences de l'indicateur de loi de puissance de 
la nouvelle fonction polynomiale de déformation de cisaillement sur la flexion de poutres à 


gradient fonctionnel sont explorées. 


2.1. Introduction 

Le concept des matériaux à gradient fonctionnel (FGM en anglais) a été mis au point 
en 1984 par des spécialistes des matériaux Japonais cherchant des matériaux super résistants à 
la chaleur destinés à être utilisés dans les navettes spatiales. Les matériaux fonctionnels (FGM) 
sont des composites microscopiquement inhomogènes souvent fabriqués à partir d'un mélange 
de métaux et de céramiques. Ils sont une classe de composites ayant une variation continue des 
propriétés du matériau à partir des surfaces supérieure et inférieure, d'une structure de poutre 
par exemple, et éliminent ainsi la concentration des contraintes à l'interface des couches 
observées dans les composites stratifiés. La combinaison de différents matériaux avec des 
propriétés physiques spécifiques permet une conception de matériau adaptée, qui étend l'espace 
de conception structurelle en mettant en œuvre une réponse multifonctionnelle avec un 
minimum de poids. La figure 2.1 présente une structure de poutre FGM dans laquelle le 


matériau passe graduellement de la partie entièrement métallique en bas à la partie supérieure 
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entièrement en céramique. La zone métallique est utilisée pour supporter les charges 
mécaniques, tandis que la partie céramique sert de protection thermique. En raison de la 
croissance des applications FGM dans les structures d'ingénierie, différentes théories de poutre 
ont été développées pour prédire la réponse mécanique des poutres fonctionnellement graduées. 
La première est la théorie des poutres classique (CBT Classical Beam Theory), connue sous le 
nom de théorie des poutres d'Euler-Bernoulli (EBBT Euler Bernoulli Beam Théory) (Euler, 
1744). C'est la plus simple des théories et elle s'applique uniquement aux poutres minces 
(élancées). Pour les poutres modérément épaisse, la CBT sous-estime la flèche et la contrainte 
de cisaillement, en ignorant l'effet de déformation de cisaillement transversale (Yang, 2008, 
Simsek, 2009, Alshorbagy, 2011). Pour cette raison elle ne s'applique pas aux poutres épaisses. 
La deuxième théorie des poutres est la théorie des poutres à déformation de cisaillement 
transversal du premier ordre (FSDBT First Shear Deformation Beam Theory). Elle est connue 
sous l’appellation de théorie des poutres de Timochenko Beam Theory (TBT Timochenko Beam 
Theory) (Timochenko, 1921, 1922). Elle a été proposée pour surmonter les limitations de la 
CBT en tenant en compte de l'effet de déformation de cisaillement transversal. Puisque la 
FSDBT viole les conditions de contraintes de cisaillement nulles sur les surfaces supérieure et 
inférieure de la poutre, un facteur de correction de cisaillement est nécessaire pour justifier la 
divergence entre l'état de contrainte réel et l'état de contrainte supposé constant (Chakraborty et 
al. 2003, Sina, 2009, Wei, 2012). Pour éviter l'utilisation du facteur de correction de cisaillement 
transversal et avoir une meilleure prédiction de la réponse de la poutre, des théories de 
déformation de cisaillement d'ordre supérieur (HSDT Higher Shear Deformation Theory) ont 
été développées et proposées par plusieurs auteurs, y compris la théorie du troisième ordre de 
Reddy (TTR Third-order Theory of Reddy). (Reddy 1984, Wang, 2000, Yesilce, 2009, 2010, 
2011), la théorie hyperbolique de Soldatos (1992), la théorie sinusoïdale de Touratier (1991), la 
théorie exponentielle de Karama et al. (2003, 2009) ; Ainsi que la formulation unifiée de Carrera 
(2003) et Carrera et al. (2011). Aydogdu (2006 & 2009) a présenté une étude dans laquelle il a 
comparé les différentes théories d'ordre supérieur, tel que la théorie de déformation de 
cisaillement parabolique (PSDT Parabolic Shear Deformation Theory), la théorie de 
déformation de cisaillement trigonométrique (TSDT Trigonometric Shear Deformation Theory) 
et la théorie de déformation de cisaillement exponentielle (ESDT Exponentielle Shear 


Deformation Theory) avec la théorie de l'élasticité tridimensionnelle. 
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Ceramic 


| Metal L s 


L 


Figure 2.1 Exemple de structure poutre en FGM 


Généralement, afin d'augmenter la précision des résultats, l’augmentation du degré du 


polynôme de l'équation cinématique est une alternative à considérer : 


U;(M, t) = U:(M €) + zw 0 (Mo, à) + z2w 2 (M, +z82w% (Mo,t) + C1) 
Avec z la coordonnée perpendiculaire au plan x-y. 

Cette technique n'est pas préférée par plusieurs chercheurs en raison de son coût élevé 
en termes de calcul. Dans ce contexte, plusieurs simplifications ont été proposées pour réduire 
le nombre de paramètres de déplacement. Une de ces simplifications consiste à simplifier les 
derniers termes de la série Taylor de la "fonction de cisaillement". La forme de la cinématique 


du mouvement proposé à travers l'épaisseur est la suivante : 
ow 
u(x, z) = up(x) — Cv + f(2)8,(%) (2.2) 


W(x,z) = Wo(x) (2.3) 
f(z) peut être considéré comme la fonction de cisaillement pour déterminer la distribution de la 
déformation et de la contrainte de cisaillement transversal à travers l'épaisseur. Selon cette 
fonction f(z), nous pouvons distinguer certains modèles, à déformation d'ordre supérieur, 
important dans la littérature, et ceux-ci sont décrits dans le tableau 2.1. 

Des théories de déformation de cisaillement d'ordre supérieur peuvent être développées 
en imposant une variation d'ordre élevé, à travers l'épaisseur de la poutre, du déplacement axial 
(Aydogdu et Taskin, 2007; Kadoli et al., 2008; Simsek, 2009; Ben-Oumrane, 2009, Li, 2010, 
Simsek, 2010, Wattanasakulpong, 2011) ou bien, à la fois, une variation d’ordre élevé des 
déplacements axial et transversal selon l’épaisseur de la poutre (1.e. en utilisant une formulation 
unifiée) (Giunta et al., 2010a, 2010b et 2011). Dans le présent chapitre nous présentons une 
théorie de poutre à déformation de cisaillement transversal d'ordre supérieur pour les poutres 
en FGM fléchies. Celle-ci est basée sur une nouvelle fonction de cisaillement polynomiale 


satisfaisant les conditions aux limites de contraintes nulles sur les surfaces libres supérieure et 
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inférieure de la poutre. Par conséquent, le facteur de correction de cisaillement n'est pas 
nécessaire. En outre, la présente théorie a de fortes similitudes avec la théorie des poutres 
classique CBT par de nombreuses expressions ; telles que les équations du mouvement, les 
conditions aux limites et les expressions résultantes des contraintes. 


Tableau 2.1. Différentes fonctions de cisaillement 


Autheurs Fonctions de Cisaillement 
Ambart 1958 j=7 I e 
mbartsumyan, f(z) = >\4" 3 
Kaczkowski, 1968; Panc, 1975; Reissner, 5z ; 47? 

1975 JO = (1-37 
47? 
Levinson, 1980; Murthy, 1981; Reddy, 1984 f(z)=7z|1-— 3p2 
; : h . [IIz 

Levy, 1877; Stein, 1986; Touratier, 1991 A rs (5) 

Soldatos, 1992 f(z) = hsin (© — Zcosh () 
| : h 2 

Karama et al., 2003, 2009 f(z) = ze-2@/h) 

: Iz 
Ferreira et al., 2005 f(z) = sin (5) 
Aydogdu, 2009 fCz) = za-26/)/na à > 0 
Mantari et al., 2011 f(z) = za 2@/R) 
Mantari, Oktem & Soares, 2012 f(z) = tan(mz) 

Il 
Mantari, Oktem & Soares, 2012 f() = sin (5 emcos(z/h) 
Z 
Mantari and Soares, 2012 f(z) = sin () ÉTECORP) 
2r 
f() = sinh1(=)-2 
Grover et al., 2013 (Go) hvr?+4 
Tr = 3 
h AT 
f() = cot”1(—)- 
Sahoo and Singh, 2013 È ) RON 
r = 0.46 


Les propriétés matérielles des poutres en FGM varient en fonction d'une distribution 


de loi de puissance de la fraction volumique des constituants. Les équations du mouvement et 
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les conditions aux limites se dérivent du principe de l'énergie potentielle totale minimale. Des 
exemples numériques sont présentés à la fin du présent chapitre pour montrer la validité et la 
précision de la présente fonction de déformation de cisaillement transversal et pour montrer les 
effets de l'indice de la loi de puissance et de la déformation de cisaillement sur la flexion des 
poutres en FGM. Des solutions numériques pour la flexion sont obtenues pour de nombreuses 


poutres sous différentes sollicitations. 


2.2. Propriétés matérielles de poutre FGM et axe neutre 

2.2.1. Propriétés matérielles effectives des poutres métal-céramique à gradients fonctionnels 
La figure 2.2 montre une poutre en FGM composée de céramique et de métal, de 

longueur L, de largeur b et d'épaisseur h. Les propriétés du matériau varient de manière continue 

et non uniforme suivant la direction z. La surface supérieure se compose seulement de 

céramique et la surface inférieure est purement en métal. Entre les deux, les fractions 

volumiques, de la céramique V. et du métal Vn, sont obtenues par une distribution de loi de 


puissance en conjonction avec la simple loi des mélanges comme suit : 


AE 2 
Va=i-v (2.5) 


Où : z est la distance de l'axe neutre ; 

p est l’indice de la loi de puissance. C’est un paramètre variable non négatif dictant le profil de 
variation du matériau à travers l'épaisseur de la poutre. La valeur de p = 0 représente une poutre 
purement en céramique ; si p est infini, la poutre est entièrement en métal, pour p = 1 la variation 
de la combinaison du métal et de la céramique est linéaire et quand la valeur de p croit, la teneur 
en métal dans la poutre FGM augmente. 

La variation de la fraction volumique de la céramique dans la direction de l'épaisseur 
des poutres FGM en fonction de différentes valeurs de p est illustrée sur la figure 2.3. Les 
propriétés effectives du matériau MPer sont évaluées par la relation : 

MPa = MP Vm(Z) + MP: V.(2) (2.6) 
où MPx et MP. représentent respectivement les propriétés matérielles des métaux et des 
céramiques. Ainsi, le module d'élasticité Eer, le coefficient de Poisson ver et le module de 


cisaillement Gerr des FGM peuvent être exprimés par les relations : 
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z 1l\P (2.7.a) 
Eerr = (Ec — En) ( " -) +En 
z 1l\P (2.7.b) 
Veff — (ve _ Vm) ( 6. à -) + Vm 
z 1l\P (2.7.c) 
Gerr = (Ge — Gm) (£ + -) + Gm 


En utilisant les expressions ci-dessus, il est possible d'obtenir une approche de la 
variation des propriétés du matériau à travers l'épaisseur de la poutre pour différents indices de 
loi de puissance. Les figures 2.4a, 2.4b et 2.4c illustrent la variation du module de Young, du 
coefficient de Poisson et du module de cisaillement d’une poutre en FGM. Il a été confirmé par 
Delale (1983), Ziou et al. (2016) et Guenfoud et al. (2016) en utilisant une méthode énergétique 
que la variation du coefficient de Poisson n'a pas beaucoup d'influence sur l'évaluation de la 
déformation. Par conséquent, nous avons adopté le même coefficient de Poisson pour les deux 


matériaux. 


Pure Alumina 


Pure Steel 


Figure 2.2 Géométrie d’une poutre FGM 


Variation de la céramique et du métal à travers l’épaisseur 


0 01 02 03 04 05 
Z/h 


Figure 2.3 Variation de la fraction volumique de la céramique suivant l’épaisseur 
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— Ceramic | 
= Ceramic 


—p-| 


Figure 2.4.c Module de cisaillement transversal G (z) 
Figure 2.4 Variation du : module de Young, coefficient de Poisson et module de cisaillement 


transversal pour une poutre FGM le long de l’épaisseur pour différentes lois de puissance 


2.2.2. Position de l’axe neutre 

La position de l'axe neutre doit être connue afin de déterminer une solution attrayante. 
Évidemment, en raison de la variation des propriétés matérielles de la poutre (le module de 
Young précisément), l'axe neutre ne se trouve plus confondu avec la ligne moyenne de la poutre, 
mais il se trouve décaler de celle-ci (ligne moyenne) à l'exception des poutres avec une variation 
symétrique du module de Young. 
Pour déterminer la position de l'axe neutre, nous considérons un nouveau système de 
coordonnées, figure 2.5, tel que le nouvel axe x est positionné sur l'axe neutre, qui sera 
déterminé comme suit. Par la suite, nous avons : 

X=X %=2+h (2.8) 

Où ho est la distance entre l'axe neutre et la ligne médiane de la poutre. Dans ce cas et de manière 


similaire au traitement habituel dans la théorie des poutres d'Euler-Bernoulli (EBBT), on peut 


= d 
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directement écrire : 


2 
Ex =-2 et ox =E(z')Ex (2.9) 


dx2 
Où w est la flèche de la poutre FGM. La position de l'axe neutre peut être déterminée en 


choisissant ho de sorte que la résultante de force axiale agissante sur la section devient nulle : 


h 
27 ho 
>r=0o/f ss dAS 0 (2.10) 
h 
En substituant l'équation (2.8) conjointement avec (2.9) dans (2.10), nous obtenons 
2 10 dw 
[ b.E(z').z'.— dz" = 0 (2.11) 
Eh, à dx 
02 


En changeant les bornes de l'intervalle, on obtient 


h 

2 dw 

2 

Alors, 
h 
d’wf fz 2 ie 
bz JF G.z de - ho [ro z |=0 (2.13) 
2 


La position de l’axe neutre se détermine par : 


h 
[2 E(z).7. dz 
hp = —2———— (2.14) 
J EC). dz 


2 


Alors, les rigidités de membrane, de couplage et de flexion deviennent : 


" ie E —En 
D. b | EG)de = Dh (E + ) (2.15.a) 
aa _h/2 m p + 1 
+h/2 bh2 p 
Dap = b | E(z)zdz = ——(E, — Eh) LT (2.15.b) 
Fr Ju D lp Dr?) 
e Er bh? p?+p+2 
Dpb = b | E(z)z’dz = (FE — En) —_î—_———_———— 15e 
li 12 (en D+3p+26 +0 Le 
En substituent les équations (2.15.a) et (2.15.b) dans l’équation (2.14) nous obtenons : 
_ se 
bo _ Fe En) +2@p+ 2 
= —————— @© — @— —— (2.16) 
h Épien 
RTE. 
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2.3. Cinématique et relations contrainte-déformation 

Considérons une poutre droite de longueur L et d'axe neutre x reliant les centres de gravité 
G de toutes les sections, avec x-z étant un plan d'inertie principal. Ici, les axes de coordonnées sont 
choisis de telle sorte que l'axe x est orienté dans la direction axiale à la ligne moyenne de la poutre 
et l'axe z positif est dirigé vers le haut et il est perpendiculaire à l'axe x, comme illustré sur la figure 
2.5. Un matériau composite FGM constitue la section transversale. Par conséquent, en général, l'axe 
moyen de la poutre ne coïncide pas avec l'axe neutre. Les charges appliquées sur la poutre sont des 
forces verticales et des moments fléchissant maintenus dans le plan x-z, comme cela est habituel 
pour les poutres planes. La flexion sur le plan y-z ne sera pas considérée ici. L'hypothèse de 
Timoshenko pour la rotation de la normale est à tenir en compte. Les déplacements axial et vertical 


d'un point quelconque A de la section de la poutre sont exprimés comme suit : 


u(x,z) = up(x) — , (x) + f(Z)b,(X) (2.18.a) 


W(x, z) = wWo(x) (2.18.b) 
Avec: 
— W,(x) est la fonction de gauchissement. 
— u(x,z) represente le déplacement longitudinal en tout point d’une section transversal. 
— wo est le déplacement transversal suivant l’axe z de la poutre. 
— up(x) = u(x,z = 0) est le déplacement longitudinal de la ligne moyenne de la poutre selon 
l’axe x. 


— f(z) est la nouvelle fonction de cisaillement adoptée dans la présente étude. 


9 97? 9 9,z\2 9 975 
ro = (5-2) = 25-20) JE 19) 
9 972 
9 
ose (©) (2.21) 


Pour le cas de la théorie classique des poutres CBT (théorie des poutres d'Euler-Bernoulli EBBT): 
f(z) = 0, tandis que pour la théorie des poutres avec déformation de cisaillement du premier ordre 
FSDBT (théorie des poutres de Timoshenko TBT): f (z) —z 


Où w% (x) est la déformation de cisaillement à la ligne moyenne de la poutre. 


x) = yx,(,z = 0) (2.22) 


La composante non nulle du tenseur de déformation est : 
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dus  d2w 0% 
ER  — 2.23 
E,(X, 2) x 22 + f(z) . (2.23) 
of 
= —, 2.24 
Yet z) Ôz LU ( ) 


G: le centre de gravité 
O: 1a position de l'axe neutre] 


Figure 2.5 La position de l'axe neutre pour une section de poutre rectangulaire 


f'(z) est nulle à la fois sur la fibre supérieure et la fibre inférieure (pour z = *h / 2 }), ce qui 
implique la nullité de la contrainte de cisaillement au niveau de ces fibres extrêmes . Il est clair que 
f''(z) est égal à zéro à z — 0. Cela implique que la contrainte de cisaillement est maximale à la 
fibre médiane. 

En introduisant l’expression de la fonction f(z) dans l’équation (2.18.a), nous obtenons 
l’expression du déplacement u en fonction de : la position du point (x, z) et la fonction du 


cisaillement. 
eo - wo. 9 97 . 
u(x,Z) = uo(x cr 87 Sn? Wy (2.25.a) 


0W 9 973 
_ er 2.25.b 
0x 8 x) L 6h2 x ( ) 


L'équation (2.25.b) a la même forme que les équations utilisées par différents auteurs tels 


que Levinson (1980), Murthy (1981) et Reddy (1984) pour les théories de déformation de 


u(x,z) = up(x) — z( 


cisaillement d'ordre supérieur (HSDT), qu’on peut la mettre sous la forme : 
u(x,z) = u,(x) + zu (x) + z2u,(x) + z$u3 (x) (2.26) 
Comme la déformation de cisaillement transversal est nulle sur les deux fibres extrêmes de 
la poutre (z = (2) ) il est donc clair que u2 (x) doit être égale à zéro. Donc, l'équation (2.26) 
+ 
devient : 


u(x,z) = uo(x) + zu, (x) + z°u:(x) (2.27) 
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Dans la présente analyse et pour des raisons de simplification, la déflection de la poutre w 
est supposée indépendante de l'épaisseur, ce qui est assez fréquent dans les théories classiques et les 
théories d'ordre supérieur des poutres (Savoia, 1996). Pour développer la relation finale de 
déformation de cisaillement transversal, nous prenons en compte en premier lieu la condition au 
limite de la contrainte de cisaillement nulle aux lignes extrêmes (fibres) de la poutre, donc la 


déformation de cisaillement y, (x, z) prend la forme 


of 9 97? 
Yxz (x, z) = dz : x = É En 7) x (2.28) 


Où w est la fonction de gauchissement (fonction de contrainte de cisaillement) de la ligne moyenne 
de la poutre. Il est à noter qu'il existe d'autres expressions satisfassent les conditions aux limites 


nécessaires. 


Figure 2.7 Gauchissement d’une section droite 


Cependant, l'expression (2.19) est l'une des expressions les plus simples, car elle se 


caractérise par une fonction polynomiale continue et régulière. 


| Tec ee ü) 
En introduisant l'angle de rotation de la section perpendiculaire à la ligne médiane 8 = = 


si Z—0; par conséquent, nous obtenons : 
Ô Wo 


2 , . ou nn : ; : 
Ainsi, la déformation normale £,(x,Z) = —. peut s’exprimée en fonction de la déflection 
X 


transversale wo et de la rotation de la section par : 


_ 0w Pwo | 9z 97 d 9z 9z3\d2w, 2.0.2) 
SOS STE le mia le 27) GS 
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, sut F9 "97 (o+2% a. 
Yxz CR Z ae 8 2hn ES (2.30.b) 


j ; Fit + : ô 
Il est claire que si Zz = 0 la déformation longitudinale devient égale à €, (x, z) = _— 
X 


2.4. Tests de validations et résultats numériques 


Considérons une poutre FGM de longueur L et de section rectangulaire. Les dimensions 
de la section sont : h=0.1m et b=h/100 = 0.001m. Deux cas sont considérés ; le rapport longueur 
sur hauteur L/h plus élevé que 100 (pour la poutre mince) et L/h inférieur à 5 (pour la poutre 
épaisse). Le module de Young suit une loi de puissance (voir équation (2.7.a)), avec Ec = 10000 
GPa et En-Ec/10=1000 GPa. Le coefficient de Poisson est égal à 0,25. La poutre est sollicitée par 
une pression uniforme égale à 1 Pa appliquée sur la fibre supérieure. Les grandeurs non- 
dimensionnelles utilisées ici sont Ec/Em = 10, G/Gm=—10. Les résultats du déplacement transversal 
(la flèche) de la poutre sont représentés sur les figures 2.8 et 2.11 pour différents exposants de loi 
de puissance p et pour différents rapports longueur sur épaisseur (L/h=—100, L/h=$ respectivement). 
Pour cette poutre en FGM, la déflection transversale augmente avec l’augmentation de l'exposant 
de loi de puissance p. Comme on le voit sur les figures 2.9 et 2.12, la répartition des contraintes 
axiales est linéaire pour le cas de la poutre en céramique pure et les valeurs des contraintes de 
traction et de compression sont égales pour le cas de la poutre isotrope (céramique complète). Mais 
pour les autres valeurs de p, la répartition des contraintes axiales n'est pas linéaire et les valeurs des 
contraintes de compression sont également supérieures aux contraintes de traction. La valeur de la 
contrainte axiale est nulle sur le plan médian (poutre isotrope et pour p=0), mais il est clairement 
visible que les valeurs des contraintes axiales ne sont pas nulles dans le plan médian de la poutre en 
FGM pour des valeurs de p non nulles. Ceci indique que le plan neutre se déplace vers le côté 
supérieur de la poutre dans le cas de la poutre en FGM. Ceci est dû à la variation du module 
d'élasticité à travers l'épaisseur de la poutre en FGM. Les figures 2.10 et 2.13 illustrent la variation 
de la contrainte de cisaillement sur l'épaisseur de la poutre pour différentes valeurs du rapport 
longueur sur épaisseur (L/h=—100, L/h=5$ respectivement) avec différentes théories de déformation 


de cisaillement d'ordre supérieur au point de coordonnée x—0. 


On peut observer que les courbes obtenues avec la nouvelle fonction de cisaillement sont 
en bon accord avec celles données par les autres théories des poutres (Kaczkowski, 1968 ; Levinson, 


1980) pour toutes les valeurs de loi de puissance p et de rapport L/h. 
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2.4.1 Premier cas : Poutre mince (rapport L/h=100) 


z/h 


Figure 2.8 Flèche W de la poutre mince en 
FGM sous charge uniforme (L/h=—100) 
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Figure 2.9 Variation de la Contrainte normale à 


travers l’épaisseur d’une poutre mince en 


FGM sous charge uniforme (L/h=—100) 
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— — kaczkowski , — — kaczkowski 


—— present present 
— — = Levinson | —— = [Levinson 
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Figure 2.10 Variation de la contrainte de cisaillement transversal à travers l’épaisseur d’une 


poutre mince en FGM sous une charge uniforme à x—0 (L/h=100) 


2.4.2 Deuxième cas : épaisse (rapport L/h=5) 


ss... p=1 


_——p2 


x p-3 
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Figure 2.11 Variation du déplacement Figure 2.12 Variation de la contrainte normale 
transversal W d’une poutre épaisse en FGM à travers l’épaisseur d’une poutre épaisse en 


sous charge uniforme (L/h=5) FGM sous charge uniforme (L/h=5) 
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Figure 2.13 Variation de la contrainte de cisaillement transversal à travers l’épaisseur pour une 


poutre épaisse en FGM sous une charge uniforme à (x—0) (L/h=5) 


2.5. Conclusion 

Nous avons développé une nouvelle fonction de cisaillement transversal pour la théorie 
des poutres à déformation de cisaillement d'ordre supérieur pour la flexion des poutres à gradient 
fonctionnel (FGM). La fonction de cisaillement transversal développée ici tient compte de la 
variation d'ordre supérieur de la contrainte de cisaillement transverse à travers l’épaisseur de la 
poutre et elle satisfait aux conditions aux limites de contrainte nulle sur les surfaces supérieure et 
inférieure de la poutre. L'utilisation d’un facteur de correction de cisaillement n'est pas nécessaire 
pour notre cas. Les solutions numériques par éléments finis obtenues avec la nouvelle fonction 
polynomiale de cisaillement transversal sont présentées aux figures 2.8 à 2.13. Cette nouvelle 
fonction a une forte convergence en comparaison avec les autres théories des poutres de 
déformation de cisaillement d'ordre supérieur (Kaczkowski, 1968, Levinson, 1980) pour divers 
exposants de la loi de puissance p et pour différentes valeurs du rapport longueur sur épaisseur. 
D'une manière générale, tous les modèles des poutres à déformation de cisaillement transversal 
d’ordre supérieur donnent des résultats différents pour la contrainte transversale de cisaillement. 


Les différentes fonctions de déformations transversales de cisaillement utilisées dans chaque 
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modèle expliquent cette différence. Les courbes de contraintes obtenues par cette nouvelle fonction 


sont régulières et continues, par conséquent le problème de délaminage ne sera pas posé. 
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ANNEXE 
Elément de poutre en FGM 


Les deux équations d’équilibre, pour une poutre en FGM à section transversal constante 
sans sollicitation, s’expriment par 


2.31 
078, : ow Sen 
E 2 By dz + GAK, (57 — 0,) = 0 


1. La Rigidité Axiale 


La matrice de rigidité axiale (membrane) [K,,] est : 


[Km] = b Î E(z)dz Es 4 


(2.32) 
2. Elément de poutre avec cisaillement en flexion 
La matrice de rigidité [K] de la poutre en flexion-cisaillement est : 
l 1 
[K] = [Ks] + [Ke] = b [ E(z)z?dz [Bsl'[Bsldx + bK, [ G(z)dz[B.]"[B.]dx (2.33) 
0 0 
Avec 
AN 
Bs] = 7 
pa = [ 
aN dy 
(Be = | “| = [2 -L 2 -1] 
d dx | 7 2L(1 + y) 
ee [ Elo? 06 060 ie 
a Ge ox 0x 
2.1. La matrice de rigidité de flexion [K/] 

12 6 —12 6 

LS L2 LS L2 
(A+2Dy+b$) —6  (2-2hy-b5) 

_ __ Dh L L2 L 
[Ks] = Gr 12 < (2.34) 
L3 L2 
(4+2Dy+63) 


É 
2.2. La matrice de Rigidité de cisaillement [K,] 


ow 0O0w / 
RE [ GG)Ky ( : o) (- 68) dv 
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Li A = À 
L 2 2 
ie L —1 L 
Dy Ds 4 2 4 
Ki] = 2, 
E 2 
L 
4 


K£et K« sont les matrices de rigidités élémentaires calculées par intégration analytique suite à 


leurs expressions polynomiales simples, elles sont données par : 


KI = Gi) 
6L = 6L 
 UY8) Le LG-à) Qu 
12 —6L 
L?(4 + y) 
2.3. Matrice de rigidité de couplage membrane-flexion 
Ke - [ro Z one 
à 0x  ôx 
La matrice de couplage s’obtient à par : 
0 0 
L. 1 
[Kepl = Gao] GG 237 
—1 1 
SR 


Par la Combinaison des différentes matrices de rigidités, nous obtenons la matrice de rigidité de 


l’élément poutre à deux nœuds : 


0 2 'APARIN EPA 0 Dab 
L L L L 
1254 64 0 1204 Où 
L3(1+dy) L2(1+dy) L3(1+dy) L2(1+dy) 
Db(4+Dy) 0 —604 Db(2-Dy) 
L(1+y) L2(1+dy) L(1+y) 
[K] = A 4 . L (2.38) 
D 0 ab 
L L 
1254 60} 
L3(1+dy) L2(1+dy) 
Db(4+Dy) 
L(1+dy) 


A Partir de l’équation d’équilibre, nous déduisons : 


Chapitre 2 Elément linéique en FGM 


ow 
0,(x) — EE + C 


Avec c est la constant d’intégration qui s’obtient à partir de l’équation d’équilibre 


028, __  GAKyc _ 12c __12/[E(z)z?dz 
ox?  BJE(z)z2d7  dyL2 Py = L2K, f G(z)dZ (2.39) 
3 
12B2 Eh 
Avec y = tD'= == 
FD 12(1—v2) 


dy est un coefficient caractérisant les déformations transversales. Il dépend des 


caractéristiques géométriques et mécaniques de la section. 
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3.1 Introduction 

Le but principal de ce chapitre et du chapitre quatre est la mise au point de la formulation 
d’un élément fini de coque permettant de prendre en compte l’ensemble des situations 
rencontrées dans le calcul des structures et ouvrages en matériau à gradient fonctionnel. La 
première étape consiste à la formulation d’un élément de coque à facette plane basé sur le 
modèle en déformation, pour la modélisation des structures composées en matériaux isotropes. 
La deuxième étape consiste à l’adaptation de cet élément fini de coque pour modéliser des 
structures composées en FGM et qui sera exposée au chapitre quatre. 

Les structures complexes en coques sont fréquemment rencontrées dans divers usages. 
Le développement d'élément fini simple et performant pour l'analyse de ces structures constitue 
une poussée majeure de la recherche scientifique en mécanique des solides. Les éléments finis 
de coque à facettes planes se dérive par superposition d'éléments de plaques avec des éléments 
membranaires de contraintes planes, pour lesquels les propriétés de flexion et de membrane 
sont découplées (pour les matériaux isotropes). Dans la littérature, il est révélé que les éléments 
de coque fournissent une bonne précision pour le déplacement et la contrainte des structures. 
Néanmoins, des problèmes sont souvent rencontrés, rendant difficile la réalisation des objectifs 
assignés. Cependant, si des éléments membranaires sont associés à des éléments de plaques en 
flexion, ils affectent la précision des résultats en raison de la pauvreté de leurs performances 
pour des problèmes de flexion dominante, ce qui nécessite une densité de maillage élevée et 
induisent des problèmes de continuité et de conformité lors de la transition vers les éléments de 
coque. Les principales contraintes observées sont souvent liées aux éléments suivants : 

- Aspects d'incompatibilité des champs de déplacement lors du couplage des éléments de 
membrane avec ceux de flexion. 

- Les deux phénomènes de blocage : de cisaillement transversal et de membrane. 

- Les problèmes numériques induits par l'absence du sixième degré de liberté dans le cas 
d'éléments coplanaires. 

- Les problèmes numériques associés à l'intégration numérique. 

De nombreux éléments finis ont été développés pour résoudre ces problèmes. Cependant, 
la plupart d'entre eux sont inefficaces pour l'analyse des configurations géométriques 
arbitraires. Les éléments iso-paramétriques sont les plus performants parmi ceux disponibles, 
en raison de leur capacité à modéliser avec succès des structures courbes. De l’analyse 
bibliographique, certains chercheurs abordent ces problèmes. Bhothikhun et Dechaumphaï 
(2014) ont développé un élément triangulaire de plaque de Kirchhoff discret (DKT) couplé avec 


l’élément de membrane CST avec maillage adaptatif. Ce dernier est utilisé pour améliorer la 
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précision de la solution et réduire l'effort de calcul. Les formulations de l’élément triangulaire 
(THS) et de l’élément quadrilatéral (QHS) basées sur le principe de la variation hybride et la 
solution homogène analytique de l'équation de la plaque mince ont été obtenues par Rezaiee- 
Pajand et Karkon (2014). Jeon et al. (2014) présentent une amélioration d'un schéma 
d'enrichissement de l'élément fini de coque triangulaire MITC à trois nœuds par des fonctions 
de recouvrement d'interpolation. Ce schéma d'amélioration augmente la précision de la solution 
sans raffinement du maillage local traditionnel. Rezaiee-Pajand et Yaghoob1 (2014) ont 
développé un nouvel élément triangulaire, nommé SST10. Sa formulation utilise les contraintes 
d'optimisation de l'insensibilité à la distorsion et à l'invariance rotationnelle. De plus, les 
équations d'équilibre sont établies en fonction de certaines conditions sur les états de 
déformations. Des éléments finis quadrilatéraux de membrane avec degrés de liberté de rotation 
ont été développés par Kugler et al. (2010), Boutagouga (2008, 2016) et Hamadi (2016). Ils 
sont basés sur un principe de variation utilisant un champ de rotation indépendant autour de la 
normale. D'autres études sur l'élément fini triangulaire ont été publiées : Sabir (1985), Barik et 
Mukhopadhyay (2002), Kim et Bathe (2009), Papanicolopulos et al. (2009), Burkardt (2010), 
Serpik (2010), Huang et al. (2010), Gileva et al. (2013) et Himeur (2008, 2015). De la même 
manière Batoz et Dhatt, en 1990, ont développé des éléments finis triangulaire de membrane 
avec trois et six nœuds appelés T3 (CST) et T6. Belarbi (2000) a également développé de 
nombreux éléments finis triangulaires, à savoir les SBT2, SBT2v, SBT3 et SBT3v. Chinosi 
(2005) a apporté plusieurs modifications aux conditions aux limites pour les plaques de 
Mindlin-Reissner. 

Récemment, Shin et Lee (2014) ont développé un élément triangulaire de coque à facette 
plane à trois nœuds basé sur la formulation en déformation déviatoire naturelle imposée pour 
une utilisation améliorée dans des structures de coques à géométrie courbe. La technique de 
l'approche en déformation est appliquée au développement de la rigidité de membrane de 
l'élément en question. 

Le but de la présente recherche, sujet du chapitre trois, est de surmonter le problème qui 
apparaît lors de la combinaison d’éléments de membrane de classe C° avec des éléments de 
plaque de classe C!. Par conséquent, la finalité de ce travail est la formulation d’un élément fini 
de coque mince à facette plane basé sur l’approche en déformation; dont le raisonnement est 
d'éviter ces difficultés d'une part et la construction d'un élément fini de coque à facette plane 
simple et compétent pour l'analyse de structures complexes, d'autre part. Pour ce faire, nous 


avons enrichi notre approche avec les concepts et le développement des techniques basées sur : 


Ka 


- l'adoption de l’approche en déformation, 
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- l'introduction d'un "quatrième nœud fictif", 

- l'élimination des degrés de liberté correspondant au "quatrième nœud fictif" par 
condensation statique, 

- l'utilisation de "l'intégration analytique" pour évaluer la matrice de rigidité. 

Dernièrement, Himeur et Guenfoud (2008, 2015) ont développé un élément fini de 
membrane triangulaire, pouvant être facilement couplé avec un élément en flexion (plaque, 
poutre et coque). Baptisé "T43_Eq" (Himeur et Guenfoud 2008, 2015) est un élément fini 
triangulaire de membrane à quatre nœuds dont un nœud positionné au centre de l’élément. Il se 
caractérise par la présence de la rotation autour de la normale, pour chaque nœud, définie par 
la dérivation des champs de déplacement (Drilling rotation). Les fonctions d’interpolation, des 
champs des déformations et par conséquent des déplacements et des contraintes, sont 
développées à partir des conditions d’équilibre (Polynômes bi-harmoniques choisis parmi les 
solutions données par (Teodorecu 1982) à partir du développement de la fonction Airy). Ces 
fonctions d'interpolation sont développées en utilisant le triangle de Pascal. L'élément, 
"T43_Eq" est un élément de membrane ayant trois nœuds positionnés au sommet du triangle et 
un quatrième nœud fictif situé au centre du triangle, tous possédant trois degrés de liberté (U et 
V et la rotation autour de la normale @,). Les degrés de liberté du quatrième nœud sont ensuite 
éliminés de la matrice de raideur élémentaire par condensation statique. La valeur principale de 
ce nœud fictif ‘esclave’ réside dans l'amélioration du champ de déplacement (raffinement P, 
c'est-à-dire augmentation du degré d'interpolation polynomiale), ce qui conduit à une meilleure 
précision dans l'approximation de la solution. Le principe variationnel utilisé est celui du 
principe du travail virtuel. L'intégration analytique pour évaluer la matrice de rigidité est 
extrêmement intéressante afin d'éviter l’insuffisance de la convergence observée lors de 
l’utilisation de l'intégration numérique dans les éléments iso-paramétriques. 

L'élément fini triangulaire de plaque en flexion possédant un quatrième nœud fictif 
(esclave) basé sur l'approche en déformation développée par Himeur et Guenfoud (2011) est 
utilisé pour construire l'élément de coque actuel. Cet élément est formulé sur la base de 
l’approche en déformation. Les fonctions d'interpolation des champs de déformations (par 
conséquent, les déplacements et les contraintes) sont développées en utilisant le triangle de 
Pascal. L'élément considéré est un élément triangulaire auquel nous avons ajouté un quatrième 
nœud fictif positionné à l'extérieur du triangle. Cette position, à l'extérieur, est donc choisie 
pour éviter la relaxation de la matrice de raideur entraînant un excès d'estimation des 
déplacements nodaux. Les degrés de liberté correspondant au quatrième nœud sont ensuite 
éliminés par la condensation statique de la matrice de rigidité au niveau élémentaire. L'intérêt 
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principal de ce nœud fictif réside donc dans l'amélioration du champ de déplacement 
(raffinement P, c'est-à-dire augmentation du degré d'interpolation polynomiale), visant une 
meilleure précision dans l'approximation de la solution. Le critère variationnel correspondant 
est le principe du travail virtuel. L'utilisation de l'intégration analytique pour l'évaluation de la 
matrice de rigidité est très intéressante pour éviter le phénomène de perte de convergence 
observé dans les éléments iso-paramétriques utilisant l'intégration numérique et qui sont très 
sensibles (leurs convergence est formée par un maillage régulier non déformé). L'hypothèse 
cinématique de la présente formulation est celle de la théorie des plaques minces (théorie de 
Kirchhoff) dont laquelle l’énergie due au cisaillement transversal est négligée. 
3.2. Formulation de l’élément de membrane (Four Node Triangular Element “T43_Eq”) 
L'élément membranaire, désigné par T43_Eq, est formulé en utilisant l'approche en 
déformation par Himeur et Guenfoud (2008). Les fonctions d'interpolation des champs de 
déformations, et donc des déplacements et des contraintes, sont développées depuis les 
conditions d'équilibre. La figure 3.1 présente cet élément. Chaque nœud a trois degrés de 
liberté : deux translations U et V et la rotation autour de la normale ("drilling rotation") 07. Par 


conséquent, les champs de déplacement ont douze constantes (4 nœuds x 3DDL-12) 


indépendantes (ar, …, a12). Les relations entre les déformations et les déplacements sont données 
comme suit : 
OV OU OV 
ri Ey = — Yo = —+— (3.1) 
ox Oy Oy x 


Figure 3.1 Elément triangulaire membranaire à quatre noeuds “T43_Eq” 
Et à trois DDL par nœud (deux translations u et v et la rotation 6,) 


La vraie rotation autour de la normale, (la "Drilling rotation") 8;, est donnée par la relation : 


MER . 
 218x ôy É 


Les conditions d'équilibre en état plan de contraintes sont données par le système suivant: 


ja) 
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Rd 
ô ô : 
se + 1-{] (3.3) 
Qu ll ’ 
oy &x|L» 


En état de contraintes planes, les relations entre les déformations et les contraintes sont 


données par les équations suivantes (loi de Hooke pour les matériaux isotropes): 


G, G 0 E, E. 1 —-v 0 . 
G,f= 12 v Î Le ë, Et E,f=—|—V E 6, (3.4) 
Le 0 0 3. Fe de 0 0 G Le 


Pour des forces de volume (nulles), les fonctions d'approximation sont développées en 
introduisant la fonction d’Airy F(x, y) qui réduit le problème de la relation (3.3) en une seule 
équation bi-harmonique. Himeur et Guenfoud (2015): 

V'F(X y) =0 (3.5) 
Avec 


co -2F&@Y , LOF, _ OF) 
x ôy’ y x: xy Exôy 


S1 nous substituons l'équation (3.6) dans les systèmes d'équation (3.4), les déformations 


(3.6) 


deviennent comme suit : 


; -1( F&Y) _, ae - 1-5 FG y). 0 mb) 


* El Gy° " à&x° E Ôy” 64 
(3.7) 
-_1{0F&y 
Lu Ox0y 


Par conséquent, la solution universelle de l'équation (3.5) se base sur les polynômes bi- 
harmoniques. Ainsi, ces fonctions sont le point de départ pour obtenir la contrainte, la 
déformation et le déplacement. Zweiling (1952) a déterminé ces polynômes bi-harmoniques. 
Pour exprimer l'approximation du champ de déformation, on utilise les douze premiers 


polynômes bi-harmoniques comme suit : 


Pour les mouvements de solide rigide (MSR), les déformations sont nulles. Donc : 


Ex =0 ni Fa (3.8) 
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BE, D EH 2.a_.V.y +2a,.x —6.a vx + 6a y — Ga, VX + Ga ;xy 


v 8 9 


BE = + 2a-.y — ZA g-wv. + Ga,.x - Ga vw. + 6a, To Ga vx 


(3.9) 


r | L _ 27. 2 
ER AE 4a,(1+v)x 4a,(+v)y Ga, (+ vx Ga (+ v)y 


Le champ résultant se caractérise par : 
- L'existence de déformations constantes garantissant la convergence et la déformation 
homogène lorsque le maillage est affiné, représentées par les constantes a4, as et a6. 
- L'existence des déformations linéaires de dilatations (& ,£,) et un état de déformation linéaire 
de distorsion (y., ) (Cisaillement dans le plan), correspondant aux paramètres a7,... ao mettant 
en dépendance les distorsions aux dilatations. 
- Les déformations dues aux modes supérieurs sont pris en considération à travers les 
paramètres a11, 12; 

L'équation universelle de compatibilité des déformations s’écrit ; 

de 82e y 


Re Mn (3.10) 
ôy? ox2  Oxôy 


L'approbation de l'équation universelle de compatibilité des déformations vérifiée par la 


fonction d’Airy s’écrit de la manière suivante : 


2 
( Fe E 62 | 1 ES y) 7 OZE(x, »] 2 1 Us y) es OTF(x, ?) 


2 
Éd y _ o {1 O2F(x, ÿ) - : O2F(x, y) l OF(x, y) N 0°F(x, y) 
6x2 6x2|El 6x? ôyZ EU &x* ôx y? 


(3.11) 


2 
. Txy __ 07 [2(1+v)| 0x y) _,20+) OTF(x, y) 
OxO0y  Ox0y E ox0y E ôx26y? 


2 2 > 
o"e, Ô Sy ê Vxy al OTE(x, y) OTF(x, y) ,9°FG, le dt : 

2 2 y E Ad ad 2 | EU ee 
Oy x ACT Ôy ox ox “y 


“E 
Finalement, l'intégration des champs de déformations fournit les champs de 


déplacement en dessous : 


D 
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_— » 262 - 2 
E.u =a a,.y+a,.x+a y 2.a V-X.y + a0.(x y=(2+v)) 3, (vx +y*) 


Il 
: 2 3 2.8 
+6a,çxy a, ,.G.vx V+Yy )+a,,.Gx Y—Y"(2+Vv)) 


gi 


x+ a? _x2(2+v))-2a 
2 


Ev=a,+a,x+ay+a 


2 
-3 ap (y +Xx )+a, .Gxy 


V.XY + 6a, X.Y (3.12) 
2 


8 
3 3 
— X (2+v))-a,,.Gvxy +x°) 


2 


Fr r Dit 2 
Eee 2.a x+2a.y+6.20.7 Ga ,x+3a (y x )+3a x”) 


3 7 
Sous forme matricielle, le champ de déplacement résultant représenté par l'équation 


(3.12) devient 


(x.y) 
vx y) ?= Lay.) (3.13) 
0, (x y) 


je 
DOME ee ne di ed 


2 2,2) 


10 -yx O0 y  -2y  (2-y/2+v) -Hx2+y 6x Cy- y) (Gx2y-y2+v) 
Dee à 3 


fix y]=l0 1 x 0 y x G2-220+v)  -2y y -3w2+x2) Gyx-xQ+v)  (Cyx-x) (3.14) 
00 1 000  -Xx 2y 6y 6x 3y2-x2) 4y2-x2) 


Ayant les quatre coordonnées nodales (x;, y;) correspondant au nœud j (j = 1, …., 4) et affectant 


la relation (3.13), le vecteur nodal de déplacements au niveau élémentaire s’obtient par : 


rx» )] | 


Le) 1 ay») Li 


=—|" (3.15) 
LESA] 
esp) 
[x.y] 
Avec {q'} =<u,,v,,8,u,,v,,67,u,,v,,8,,u,, V8, > et — ee est la matrice 
[x 4.7.) 


élémentaire des coordonnées nodales qu’on la note par [A]. 
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Himeur et Guenfoud (2008, 2008, 2015) détaillent le développement de la matrice [A] 
pour l'élément fini "T43_Eq". Enfin, nous obtenons les fonctions d'interpolation des 
déplacements à partir de l'équation (3.15) par l'évaluation des paramètres ‘ai’ en fonction des 
coordonnées nodales. 

{a }=[AT {a} (3.16) 

En substituant l'équation (3.16) dans l'équation (3.14), nous obtenons la relation entre le 
vecteur des déplacements et les déplacements nodaux par la relation suivante : 

(x.y) s 
vG y f= Lx ar tfe) (3.17) 
0,(% y) 


sl x.y) ITA] 7 =[N] de taille 3x12 représente la matrice des fonctions de forme. 


Dans le cas de l’état plan des contraintes, le tenseur de déformation de Green linéarisé 


représentant les relations déformations-déplacements est donné par : 


eu 
ul SES 
u(x, y 
£ =| 0 —| 3.18 
M Oy ei 6518) 
1e 0 à 
7, er. 2e 


En remplaçant u et v par leurs valeurs de l'équation (3.17), les relations déformations- 


déplacements de l’équation (3.18) prennent la forme suivante : 


Ex 0 0 0 1 O0 0 - 2vv 2x -6vv 6y - Evvx 6xy 
Bis, |=l0 0 0 0 1 0 2 _2w 6x -6w  6xy 6x fe} (3.19) 
Yy] L0 0 0 0 0 2 a+ -414+n 0 0 _6x2(1+v) -6y2(1+v) 
Par conséquent, la matrice de déformation prend la forme : 
0 0 0 1 O0 0 -2vv 2x -6vv  6y - GEvVx 6xy 
(o]- 0 010: 0 0 te : Ge. EON (6 -Gwx | (3.20) 
0 O0 0 0 O 2 -4x(1+v) -4y(1+v) O0 O0  -6x2(1+v) -6y2(1+v) 


3.3. L’élément de plaque 

3.3.1 Les équations fondamentales de la théorie des plaques minces (Théorie de Kirchhoff) 

3.3.1.1 Les équations de la cinématique des plaques minces et les conditions de compatibilité 
La figure 3.2 présente les rotations autour des deux axes x et y désignées par 6, et 6, 

respectivement et les pentes (gradients) dans les deux directions définies par les variables f, et 


By; avec: 


RE 
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Bx=0) By =-0% (3.21) 


Deformed 
midsurface 


Original 
midsurface 
x 


Section z = 0 


Figure 3.2 Déformation d'une plaque en flexion (Théorie de Kirchhoff) 


L'hypothèse de la section transversale plane avant transformation reste plane après 
transformation, engageant ainsi une variation linéaire du déplacement à travers l'épaisseur de 
la plaque, se traduit par la relation (3.22). 
u(x,y,2)= 20; (a, y)= 10, (y) v@7,27)= 170, (a y)=-10, (x, y) WG y,2)=w0 y) (3.22) 


Les trois expressions de l’équation (3.22) permettent de relier les déplacements 
membranaires (u, v) et le déplacement transversal de flexion (w) en référence à l'hypothèse de 


Kirchhoff. Le champ de déplacement final nécessaire pour décrire le comportement de la plaque 


devient : 
U = —Z Ow V = —Z Ow 2 
x Ce (3.23) 
Les relations entre les rotations, les courbures et les pentes sont données par: 
Ow Ow 
“6x =f6y 7 © CAN PR (3.24) 


Le tenseur de déformation infinitésimal s’obtient par dérivation du champ de 


déplacements, ce qui nous donne : 


p M 0x Ow , _Ov_ By __ 87w 
7 Ox &x2 y Ôy ©y oy2 
; (3.25) 
Ou . Ov 0Bx By O=w 
Vxy — = Z( ) = —22Z 
Oy x Oy x OxOy 
AVEC Yxz = Yyz = 0, Suite aux hypothèse de Kirchhoff. 
Les relations entre les moments de flexion et les courbures s’expriment par : 
OBx __07w By __ o7w OBx , By o2w 
K y; Ky= — K xy =( + )= 2 (3.26) 
ox a2 Oy oy2 Oy x ©xôy 


Afin d’assurer l’unicité de la solution, Saint-Venant en (1854) a établi les conditions de 


compatibilité. Dans le cas des plaques minces elles s’écrivent comme suit : 


3 
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2 2 2 2 
O7Ky ° K _Ô Kxy 6? yye É Vyz di Kxy _9%K% 
ôy2 dx? Ex0y OxÔy dx? êx ©y 
2 2 | 
D'Yye _Oyyz , OKay _ OK y 
Ox0y ôy2 ©y üx 


3.3.1.2 Lois Consfitutive 
La loi de comportement, dans le cas de l’état plan des contraintes et pour un matériau 


isotrope, s’écrit : 


Ox . 1 v O0 ||e, 
Ti 0 0 2. Vxy 


_27w 

M, 1 v Oo |[K, M, 1 v O0 à? 

Eh Eh o2w 

My f=—— |" 1 Ni K, Myf=——|v 1 0 É-— (3.29) 

Dev") = 12(1-v2) 1-v|| y 

M xy 0 0 —|K,, M xy Due 2, 
=0 

Ox0y 


3.3.1.3 Les équations d’équilibre 
L'équilibre entre les actions internes et externes d'un élément de plaque, figure 3.3, de 


dimension dxxdy s’écrit : 


0Q 
qdxdy +(Q, + _ dx )dy + Q, _ dy)dx -Q, dy Q,dx =0 (3.30) 


Où Q, et Q, sont les efforts tranchants (les forces de cisaillement transversal) dans les sections, 


perpendiculaires aux axes x et y respectivement; l'expression (3.30) se simplifie pour devenir : 


6Qx 20 
ee (3.31) 
L'équilibre rotationnel des moments autour des axes x et y fournit : 
M y M y My My 
= ER, 32 
eme 0 a à 2 


En substituant les valeurs de Qx et Qy exprimés par les équations (3.32) dans l'équation 
(3.31) et en utilisant la loi de comportement de flexion (3.29), la condition d'équilibre se traduise 
en fonction du déplacement transversal "w" par la relation suivante : 
04w 04w " 04w q _ 
ôx* x 26y2 ôy4 D 


0 (3.33) 
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Figure 3.3 Equilibre des forces et des moments sur un élément de plaque 


3.3.2 L'élément fini de plaque mince développé par Himeur et Guenfoud 
3.3.2.1 Les fonctions de forme 
Pour les mouvements de solide rigide (MSR), les courbures sont égales à zéro 
Kx =0 Ky = 0 Kxy = 0 (3.34) 
En changeant les courbures avec leurs valeurs données par les équations (3.34) dans les 
équations (3.26) et ensuite par intégration, nous obtenons les champs de déplacement 


correspondant aux mouvements de corps rigides, qui sont donnés par : 


W=a,-a,X-a.y B,=a, Ponts (3.35) 


Avec à et az, les paramètres liés aux rotations 6, et 0, du corps rigide autour des axes 
LL PS 


y" et "x" respectivement et a1 (le déplacement transversal) du corps rigide selon la normale 


(axe "z") lié à la translation transversale. 


L'élément possède quatre nœuds (Figure 3.4) : trois en tête du triangle auxquelles nous 
rajoutons un quatrième nœud virtuel (auxiliaire). Tous les nœuds ont trois degrés de liberté. Les 
champs des déplacements, exprimés par l'utilisation du modèle en déformation, nécessitent 
donc douze paramètres indépendants (ar, …, a12). Les trois premiers (ar, &, a3) sont utilisés dans 
(3.35) pour représenter les mouvements du solide rigide (MSR). Les neuf restants (a4, …, a12) 
sont utilisés pour caractériser l'état de flexion pure. Ils interviennent dans les fonctions 


d'interpolation de déformation pour satisfaire les relations (3.27) de compatibilité cinématique 


REA 
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pour l'élasticité plane et les équations d’équilibre. Par conséquent, les champs de courbure de 


flexion pour les modes supérieurs dérivant du triangle de Pascal sont : 


K, ay tas X+ae -V+ATX Ë, =ag +agX+a[py +4 x 


_ 2 
=) +246X+ ax +249.) +4 


; (3.36) 
11Y 


y 


+ Un quatrième nœud fictif 


Figure 3.4 Elément triangulaire de plaque mince avec trois degrés de liberté par nœud 


En substituant, dans les relations (3.26), les courbures avec leurs valeurs données par 


les équations (3.36) et en intégrant les équations résultantes, on obtient le champ de 
déplacements du mode de flexion pure : 


W = x? x? x2. x. y2 x.y2 y? x. y9 x.y 
a 4. ÿ “As. G “A6: à G Bag 9 5 40 G CITE 6 GT 
2 2 3 
— J X J J : 1 J y J y 
Bx=ayx+as. : HaÇ-x) +47. 5 ao: : +ay : *a25 (3.37) 
2 3 2 2 
— X J X J J J y L X.y J X 
By a6- 3 a. ci ag-Y+49XY+4j0- 5° a] . 25 


En additionnant les relations (3.35) et (3.37), on obtient le champ de déplacement final 


(correspondant aux mouvements du solide rigide et à la déformation pure en flexion) : 


2 3 2 3 2 2 3 
: X x.y x”: x.y 
W dj—dpX— A3) a 4 2 “As: 6 “A6. 3 a. 6 dg- n do- 5 dj0: 6 
x.y? x.y 
1 AS 
) __. (3.38) 
* } Y: y J y 
Bx=a; ayX+as ax.) +a. +4) - +ajj- 6 +-aj) - 
2 3 2 2 
By=a ne +ao.Y+ägx.y+a _. +, - 
Y73 ‘6 T6 8 9 10 


+ d1 1. rad ; 
2 11° 9 12° 


Le champ de déplacement donné par les équations (3.38) s’écrit sous la forme 
matricielle comme suit : 
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WG, y) 
p œy=[ix ph.) (3.39) 
B &y 
y 
AVE : 
km =<a à ,4 ,4 ,4_,4 ,4_,4 ,4 ,4 ,4 ,4 > 
1 123 4 5 6 7 8 9 10 11 12 
Et 
1 e y x2 x3 XT.y X7.Y y2 x.y2 y? x.y9 X.Yy 
2 : 2 5° 2 z 6 ? 2 
: x X°.y y y y 
F&y]=)0 1 0 x — x) ; GO = - 6, > (3.40) 
2 3 
A SR SO : 
2 6 2 2 2 


Le vecteur élémentaire des déplacements nodaux correspondant aux nœuds j (j = 1 … 4), 


s’obtient en appliquant la relation (3.40) après connaissance des coordonnées nodales (x;, y;) : 


GX) 
f(x, ,Y) 
We) 27218.) (41) 
(EPP 

EX 4: V4) 
F7] 
T f(x, y) 
Avec, gel = Grp Po ya Ps y g Bray et WALES 
3° 73 
EX 4 V4) 


L’expression détaillée de la matrice des coordonnées nodales se trouve dans Himeur et 
Guenfoud (2011) en Annexe. 
À partir de la relation (3.41), nous recueillons les valeurs des paramètres "ai" par le 


système suivant: 
h,J=rar fe) (3.42) 


En remplaçant les paramètres exprimés par la relation (3.42) dans le système d'équation 
(3.39), on obtient la relation reliant le vecteur de déplacement au vecteur des déplacements 


nodaux : 


REA 
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WX, y) 
= = Q 

B, y)? = fx, y) HAI ÿ | (3.43) 

BC y) 

y 


Donnant ainsi la matrice des fonctions d'interpolation Ni pour le présent élément de 
plaque mince de Kirchhoff. 
En remplaçant dans les équations (3.26), w (x, y) par les valeurs de l'équation (3.39), La 


relation déformation-déplacement prend la forme développée suivante 


S 0 O0 O 1 x y xy 0 0 0 0 0 

K}= Ke (=0 000.0: 0 0 1€ 5 of. (3.44) 
K 0 0 0 0 0 2x x? O0 2y O y? 1 
xy 


Finalement, la matrice de déformation [Q(x, y)] de l’élément de plaque en question 


s'obtient par : ([K] = [Q (x, y)I {ai}): 


0 O0 O0 1 x y xy 0 0 0 0 0 
[QG »]=10 0 0 0 0 0 0 1 x y xy 0 (3.45) 
0 0 O0 O0 O0 2x x2 0 2y 0 y? 1 


3.4. Les matrices de rigidité élémentaire pour les deux éléments (membrane et plaque en 
flexion) 


Le travail virtuel interne élémentaire discrétisé s’exprime par : 
d T 
(ow, ) = {6e}! {o}dv (3.46) 


Sachant que : 
{E}= EN Ha°}= [QG »HTAI ge} (3.47) 

Ft t}=[Dlt} (3.48) 
De plus, la substitution dans l'expression (3.46) de {} et {o} par leurs valeurs données 


respectivement dans les équations (3.47) et (3.48) produit : 
€ 
(ow, ) = fe} Ja {ec nf Dec tal “ie (3.49) 
V 


Par conséquent, on obtient la matrice de rigidité élémentaire dérivée de l'expression 
(3.49), à la fois pour les éléments de membrane et de flexion en utilisant l’intégration sur le 


volume de l’élément : 


Re 
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[K°]= ra doc pITIDIQx y ITAT lav (3.50) 


L'expression (3.50) peut s’écrire comme suit : 


[R°1=[a TT | cp} DifQc pjavrart 214 Tir AT (3.51) 


L'évaluation de l'expression de [Ko] se détermine par intégration analytique Belarbi 


(1999), Hamadi (2006), Himeur (2008), des différents composants du produit matriciel 


résultant [Q(x, y DI.[Q(x, |. Dont leurs expressions prennent la forme " Hg > Cr yP ". Pour 


plus de détails voir l’annexe IV. 

En conclusion, la matrice de rigidité élémentaire, à considérer lors de l'assemblage et de 
la construction de la matrice de rigidité globale de la membrane et de la plaque en flexion, est 
acquise après condensation de la matrice [Ke]. La condensation statique concerne les degrés de 
liberté du quatrième nœud fictif, Himeur et Guenfoud (2008, 2008, 2011, 2014 et 2015) et 
Guenfoud (2018). 

3.5. L’élément de coque 

Le présent élément est un élément de coque mince à facette plane, obtenu en superposant 
l’élément de membrane ‘T43_Eq” Himeur and Guenfoud (2008) avec l'élément fini de plaque 
mince Himeur and Guenfoud (2011). Pour le cas de matériau isotrope, nous obtenons sa matrice 
de rigidité élémentaire en ajoutant la matrice de raideur de l'élément de membrane à celle de 
l'élément de flexion sans effet de couplage Guenfoud (2018). Nous décrivons l'approche utilisée 
avec ses fondements comme suit : 

* Nous approximons la géométrie réelle de la structure par des éléments à facettes planes 
(Figure 3.5), donc nous négligeons les courbures sur l'élément et par conséquent nous évitons 


le blocage de la membrane. 


Figure 3.5 Elément de plaque mince dans les deux référentiels (global et local) 
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* Nous utilisons un élément de membrane couplé à celui d'un élément de plaque en flexion 


(Figure 3.6). 


3 
1 1 
Elément Membranaire Elément de plaque mince Coque à facette plane 
9 DDL 9 DDL 18 DDL 


Wu, v, & u, v, @&, u, v, be, 0, 0, w 6, 6, w 6, 8. 18 DDL 


T43_Eq Elément de plaque mince Elément 


de coque 


Figure 3.6 Mode de construction de l’élément de coque 


* L'élément de coque peut avoir n'importe quelle orientation dans le système de coordonnées 
global XYZ (Figure 3.5). 
* Le passage des coordonnées locales aux coordonnées globales, pour l’évaluation de la matrice 


de rigidité élémentaire dans le système global, se fait à travers la matrice de rotation [Ro] comme 


suit : (3.52) 
1] AIM (3.52) 


* Nous réarrangeons les termes de la matrice de rigidité locale (18x18) avant assemblage 


comme sur la Figure 3.7. 


Avec: 
Terme de rigidité de membrane | 
Terme de rigidité de flexion EL] 


Nous levons la difficulté associée à l’absence de la rigidité autour de 8;, dans la formulation de 
l'élément membranaire, en introduisant la vraie rotation autour de la normale «Driling rotation» 


dans la construction de la matrice de rigidité élémentaire. 


Eu 
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Nœud 1 Nœud 2 Nœud 3 


U [Vi Wilôxi Ov: Oz | U | V2 | W2|0x Ov2 02 | U3 | V3|W3) 0x3 0y3 073 
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Vi 
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Oz 
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Va 
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Ov2 
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0x3 
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Figure 3.7 Structure type de la matrice de rigidité élémentaire dans le repère local [Ke] 


3.6. Validation de l’élément coque — cas de matériau isotrope- 

Nous examinons un certain nombre de problèmes de référence, couramment utilisés par les 
chercheurs, pour comparer le présent élément avec d'autres se trouvant dans la littérature afin 
de mettre en valeur sa précision relative et sa convergence. Nous récapitulons comme suit les 


éléments finis utilisés dans les comparaisons : 


EE 
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ANST3 


ANST6 


ECBI 


ECB2 


DKTM 


DSTCOQ 
DSTM 


DKTI2 
DKTI8 
DLRI8 
ACM_RSBES 
ACM_SBQ4 
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Description de l’élément 

Elément fini de coque triangulaire à trois nœuds basé sur les 
déformations covariantes imposées Guenfoud (1990). 

Elément fini de coque triangulaire à six nœuds basé sur les 
déformations covariantes imposées Guenfoud (1990). 

Elément fini de coque triangulaire basé sur la théorie des coques 
profonde Geoffroy (1983) 

Elément fini de coque triangulaire basé sur la théorie des coques 
surbaissé Geoffroy (1983) 

Elément fini surbaissé de coque mince à trois nœuds basé sur la 
théorie de coque surbaissée de Marguerre El-Khaldi (1987). 
Elément fini de coque à facette plane Guenfoud (1996, 2000). 
Elément fini de coque épaisse à trois nœuds basé sur la théorie de 
de coque surbaissée de Marguerre Guenfoud (1996, 2000). 
(CST+DKT6) Batoz (1990). 

(CST+DKT9) Batoz (1990). 

Carpenter (1986). 

ACM_RSBES : Hamadi (2016). 

ACM_SBQ : Belarbi (2000). 


3.6.1 Analyse spectrale de la matrice de rigidité élémentaire 


L'évaluation des valeurs propres de la matrice de rigidité du présent élément est réalisée 


sur différentes formes d'éléments. En effectuant l’analyse spectrale de la matrice de rigidité du 
présent élément, nous révélons six valeurs propres nulles associées aux six modes de 
déformations du solide rigide requis. Aucun mode d'énergie parasite nulle n'est présent, comme 
on aura s'y attendre d'une matrice de rigidité intégrée analytiquement et exactement. Par 
conséquent, la matrice de rigidité a un rang correct et l'élément est considéré comme cohérent 


sur le plan cinématique. 


3.6.2 Patch Test 
3.6.2.1 Patch test standard de Kirchhoff 

Nous considérons le patch test standard de plaque de Kirchhoff décrit comme sur la 
figure 3.8. Une plaque mince rectangulaire, soumise à un état de flexion, de sorte que tous les 


résultats des trois moments soient constants sur tout le domaine de la plaque. La plaque est 


LA 
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discrétisé en quatre éléments géométriquement distincts, produit des valeurs exactes des trois 


moments et elles sont constantes sur tous les points de l'élément. 


E=1000 h=1 


v=0.3 
Boundary restraints 
W=( at nodes 1, 2 and 4 


Figure 3.8 Patch-test standard de Kirchhoff 
3.6.2.2 Flexion inextensible 
Afin d'étudier le problème de blocage (le verrouillage) de membrane, nous considérons 


une plaque rectangulaire initialement courbe soumise à une flexion constante. Dans le premier 


4Hx 


— (a = =), tandis que pour le second 


cas, la flèche initiale est donnée par l’expression W, = - 


cas elle est exprimé par Wo=Hxy/BL. Nous analysons le test pour une épaisseur H= 0.1 et 1. La 
solution analytique, selon la théorie des coques surbaissée de Marguerre, est donnée par 
Guenfoud (1993). Dans tous les cas, l'erreur sur l’évaluation de l’énergie de déformation est 
inférieure à 1,1%. La contrainte maximale de membrane est égale à 2,9% de la contrainte de 
flexion pour H=1 et inférieure à 1,1% de la contrainte de flexion pour H=0, 1. Nous restreignons 
que l'élément ne se verrouille pas. Si nous n’avons pas fait une attention particulière au 


verrouillage de la membrane, les résultats seraient très modestes. 


3.6.3 Cylindre pincé à bords libres 

Le cylindre pincé à bords libres est un problème classique qui a été largement utilisé 
pour vérifier la capacité des éléments de coque à représenter une déformation inextensionnelle 
dans le cas de flexion. C’est un exemple à utilisation courante dans la littérature. En effet, depuis 
1967, il a toujours servi aux chercheurs, Thomas (1975), Zienkiewicz (1977), Batoz (1977), 
Geoffroy (1983), Carpenter (1986), Guenfoud (1990, 1996), Guenfoud (2018), pour évaluer la 
performance de la convergence des éléments finis de coque d’une part et de vérifier la 
représentation du mouvement du corps rigide d’autre part. Le test concerne l'analyse d’un 


cylindre pincé isotrope à bords libres se trouvant à la figure 3.9. Les extrémités libres du 


Re 
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cylindre conduit à une déformation de flexion pure sans extension lorsque h/R se rapproche de 
zéro. Pour le présent exemple, nous avons h/R=0,0031 et h/L=0,0015, indiquent que le cylindre 
est très mince. On représente sur la figure 3.9 la géométrie et les données mécaniques du modèle 
à analyser. Le chargement concerne l'application de deux forces diamétralement opposées (P) 
et (-P) sur la section médiane du cylindre. 

Tableau 3.1 Cylindre pincée à bords libres 


— Caractéristiques géométriques et mécaniques- 


Longueur du cylindre ( L) 10.35 
Rayon du cylindre(R) 4.953 
Module de Young (E ) 10.50x106 
Coefficient de Poisson (v) 0.3125 
Epaisseur (h) 0.01548 
Force concentrée (P) 0.10 


Du fait de la double symétrie, nous modélisons seulement un huitième ABCD du 
cylindre en imposant les conditions de symétrie le long des côtés "AB", "CB" et "DC". Les 
conditions aux limites considérées ici sont les suivantes : 

Tableau 3.2. Cylindre pincée à bords libres 
— Conditions aux limites et condition de symétrie- 
W=6;=0,=0 Along AB 
V=60x;=0;=0 Along CB 
U=8,=0,;=0 Along DC 


Triangular 
elements 


Figure 3.9 Cylindre pincé à bords libres 
Le maillage considéré est composé d'une maille unique (double élément triangulaire) 
suivant la direction longitudinale, puisque la coque se déforme sensiblement de la même 
manière dans chaque section transversale, et plusieurs mailles sont considérées suivant la 
direction circonférentielle. Les rapports L/h=668,6 et R/h=320 indiquent que la coque est très 


mince. Timochenko et Woinowsky-Krieger (1959) donnent une solution analytique pour ce cas 


EI 
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limite (coque mince profonde), qui est Wc=-0,0244. Nous comparons nos résultats avec ceux 
donnés par Guenfoud (1990, 1996). Nous notons que pour le présent exemple, le mode de 


flexion est dominant. En effet, au point B, par exemple, nous avons pour z=h/2 que : 


N 12 h N 12 À 
Gx(2) = + M,..2 GZ =D) = + MES 


D'oùs Æ=°} 
Nx 6 

Le paramètre considéré ici pour évaluer les performances du présent élément fini de 
coque pour l'exemple actuel, sont le déplacement suivant l'axe z du point € (Le déplacement 
transversal Woc) et le déplacement radial au point B (déplacement UB). Nous présentons aux 
figures 3.10 et 3.11 les résultats de notre élément et nous concluons ce qui suit : 
- Une convergence rapide du présent élément, reflétant ainsi une représentation correcte de l'état 
des mouvements du corps rigide, 
- Notre élément converge vers la solution donnée par la théorie des coques profondes Wc=- 
0.0244 (cas du cylindre mince). 
- Avec le maillage sélectionné, nous observons une convergence monotone ascendante du 
présent élément vers la solution de référence. 
- Les résultats numériques obtenus pour le cylindre pincé à bords libres indiquent que l'élément 


développé présente une bonne précision. 


Pinched cylinder with free edges 
Diplacement W at point € 
0,025 
0,0245 
0,024 
0,0235 


0,023 —+— ANST3 
Z 0,0225 —— ANST6 

0,022 Exact 

0,0215 

0.021 —€— Present 

0,0205 ——DSTCOQ 


0.02 
40 90 140 190 
Number of DOF 


Figure 3.10 Cylindre pincée à bords libres — Déplacement W au nœud C - 
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Pinched cylinder with free edges 
Diplacement U at point B 


0,023 
0,0225 


0,022 Bu" 


0,0215 


4 


2 0,021 
0,0205 
0,02 


0,0195 


: 


0,019 
40 90 140 190 
Number of DOF 


Figure 3.11 Cylindre pincée à bords libres — Déplacement U au nœud B 


3.6.4. Cylindre Pincé Infiniment long 

Le cylindre circulaire pincé infiniment long, exposé sur la figure 3.12, est soumis à deux 
lignes de charge transversales uniformes diamétralement opposées. Nous cherchons à résoudre 
ce problème de manière à démontrer l'influence des relations entre déformations et 
déplacements d'une part et à suivre le comportement de notre élément tout en modélisant des 
structures courbes d'autre part. Nous présentons les données géométriques et mécaniques de 
cette coque sur la figure 3.12. De la même manière que dans Batoz (1977), Geoffroy (1983), 
Guenfoud H. (2017), Guenfoud H. (2018) et Guenfoud H. (2018) et sur la base de la symétrie, 
nous ne modélisons qu'un quart EAFDBC d'une bande de la coque. Les conditions aux limites 
imposées permettent de se libérer de l'influence de la direction Y sur la déformation de 
l'ensemble. Des forces ponctuelles substituent la charge uniforme. 

Tableau 3.3. Cylindre Pincé Infiniment long 


— Caractéristiques géométriques et mécaniques- 


Rayon du cylinder (R ) 10 
Module de Young (E) 10° 
Coefficient de Poisson (v) 0 
Epaisseur (h) Il 


Tableau 3.4. Cylindre Pincé Infiniment long 


— Conditions aux limites - 


V=6;=0 En tous noeuds 
U=0,=0 Le long EF et DC 


[#) 
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Figure 3.12 Cylindre Pincé Infiniment Long 


Nous dévoilons les résultats numériques normalisés des déplacements verticaux pour 


les nœuds A et B en fonction du nombre total de degrés de liberté sur les Figures 3.13 et 3.14. 


Nous comparons les résultats de notre élément à deux solutions exactes. La première, donnée 


par Donnel (1933), elle est basée sur la théorie des coques surbaissées et la seconde basée sur 


la théorie des coques profondes donnée par Koïter (1960, 1973) et Sander (1959). De même 


que dans Guenfoud (1990) et Geoffroy (1983), Il est à noter que : 


- Une convergence monotone rapide vers la solution exacte est obtenue par notre élément. 


- On remarque l'existence d'un saut entre les deux premiers maillages et le troisième. Ceci est 


dû au fait que la coque a un comportement de coque peu profonde pour les deux premiers 


maillages alors qu'elle montre un comportement de coque profonde lorsque le maillage est 


raffiné. 
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Figure 3.13 Cylindre Pincé Infiniment long 
Déplacement Normal au nœud A 
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Infinitely long pinched cylinder 
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Figure 3.14 Cylindre Pincé Infiniment long 
Déplacement Normal au nœud B 


3.6.5 Panneau Cylindrique soumis à son poids propre (Scordelis-Lo roof, coque mince peu 
profonde 1964) 

Le prochain test, considéré ici et utilisé de façon récurrente pour étudier la performance 
des éléments de coque, est le panneau (toit) de Scordelis-Lo (1964). C'est certainement le 
problème le plus utilisé pour comparer les différents éléments de coque proposés dans la 
littérature. C'est un panneau cylindrique circulaire où les deux bords courbes sont simplement 
supportés par deux diaphragmes rigides le long de leurs côtés plans et les deux autres bords 
sont libres. Le panneau est soumis à son propre poids seulement. Nous donnons à la figure 3.15 


et au tableau 3.5 les caractéristiques géométriques et mécaniques du panneau. 


Tableau 3.5. Panneau Cylindrique soumis à son poids propre 


— Caractéristiques géométriques et mécaniques- 


Longueur du Panneau ( L) L=6m 


Rayon moyen de la partie courbe (R) R=3m 


Epaisseur du panneau (h) h=0,03m 
Module de Young (E) 30000N/m° 
Coefficient de Poisson (v) 0,3 


Poids volumique du panneau (5) Ve=-0,2083N/m° 
Angle p—40° p =40° 
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Tableau 3.6. Panneau Cylindrique soumis à son poids propre 


— Conditions aux limites et conditions de symétries - 


Conditions aux limites 
U=W=0,=0 Le long de AD 
Conditions de  symétries 
V=6;=0,=0 Le long de CB 
U=0,-0,=0 Le long de CD 


Tableau 3.7. Panneau Cylindrique soumis à son poids propre 


— Valeurs de comparaisons - 


Valeurs de références _ Théorie des coques profondes 
WB=-0,00361 m Wc=0,000541m 
Solutions Analytique Théorie des coques surbaissées 
WB=-0,003703 m Wc=0,000525m 
Ug=-0,001965 m Va=- 0,001413m 


Rigid disphragm\ | 


Figure 3.15 Panneau Cylindrique soumis à son poids propre 


Scordelis-Lo roof, Coque mince surbaissée 


L’actuel panneau est mince et peu profond (R/h = 100, L/h = 200). La symétrie du 
problème nous permet d'étudier seulement le quart ABCD de la structure avec des maillages 
réguliers (N = 2, 4, 6 … éléments) aux côtés de AB et AD. Nous négligeons les déformations 
transversales dues au cisaillement transversal. Les déformations membranaires (dues aux 
efforts normaux) sont plus importantes en comparaison avec les déformations de flexion (dues 


aux moments fléchissant). L'influence de la membrane sur le comportement global de la coque 
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est très dominante. En effet, la comparaison des deux effets rapportés (membrane et flexion) au 


hNy 


= 4, 
6Mx 


point (B) donne : 


Les résultats obtenus par le présent élément de coque à facettes planes, développé dans 
cette thèse, pour le déplacement vertical au point médian B du bord libre et au centre C du 
panneau (toit) sont montrés aux Figures 3.16 et 3.17 respectivement. Nous comparons nos 
résultats à deux solutions de référence. La première analytique repose sur la théorie des coques 
profondes due à Scordellis (1964) et la seconde basée sur la théorie des coques surbaïissées due 
à Forsberg (1970) d'une part et aux résultats numériques de la littérature d'autre part. Nous 
constatons une convergence rapide des résultats du présent élément vers la solution des coques 
peu profondes par rapport à d'autres éléments. Cette convergence remarquable due à la richesse 
de l'élément en membrane (interpolation cubique). Les courbes de convergence des Figures 
3.16 et 3.17 montrent la bonne contribution de l'approche en déformation. Les résultats 
généraux approuvent la bonne convergence de l'élément de coque formulé dans le présent 


travail Guenfoud H. (20138) et Guenfoud H. et al. (2018). 
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Figure 3.16 Panneau Cylindrique soumis à son poids propre 


Déplacement Vertical au nœud B 
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Scordelis-Lo roof 
Vertical Displacement at point C 
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Figure 3.17 Panneau Cylindrique soumis à son poids propre 
Déplacement Vertical au nœud C 
3.6.6 Coque sphérique à base carrée 
Le problème présenté sur la figure 3.18 est la coque sphérique à base carrée articulée 
sur son pourtour. Cette coque se caractérise par son rayon R, son ouverture a et son épaisseur 
h. Elle est soumise à une charge concentrée (P) au centre de sa crête. La figure 3.18 présente 
les données géométriques et les caractéristiques mécaniques de cette coque sphérique à base 


carrée. Tous les bords de la coque sont articulés. 


Tableau 3.8. Coque sphérique à base carrée 


— Caractéristiques géométriques et mécaniques- 


Longueur de la base carrée (2a) 2a=1569,8 mm 


Rayon de la sphère (R) R=2540 mm 
Epaisseur (h) h=99,45mm 
Module de Young (E ) 68,95 MPa 
Coefficient de Poisson (v) 0,3 

Force concentrée (P) P=IN 


Es 
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Figure 3.18 Coque sphérique à base carrée 


Données géométriques et maillage 

Les rapports 2a/h=15,78 et R/h=25,54 indiquent que la coque est modérément mince. 
Elle est également assez surbaissé car le rapport de sa profondeur à sa longueur de côté est égal 
à 2a/H= 12,65. La symétrie du problème nous permet d'étudier seulement le quart ABDC de la 
coque. Ce problème a été traité analytiquement par Leicester (1968) et numériquement par 
plusieurs chercheurs tels que Dhatt (1970), Gallagher (1975), Tahiani et Lachance (1975), 
Fezans (1981), Guenfoud (1990, 1996) et Guenfoud H. (2017, 2018). Les résultats obtenus pour 
les déplacements verticaux au point de chargement (C) sont montrés sur la figure 3.19. La 
solution de référence est tirée de Fezans (1981). On observe une remarquable convergence de 
notre élément pour cet exemple, ce qui confirme l'utilité de notre élément pour la modélisation 
des structures à double courbure surbaissées. Aïnsi, il converge de manière monotone vers la 


solution des coques peu profondes. 


Square base spherical shell 
Displacement at point € 


Deflection at point C 


20 120 220 320 420 520 
Number of DOF 
Figure 3.19 Coque sphérique à base carrée 
Déplacement vertical au centre € de la coque 
3.6.7 Coque Parabolique Hyperboloïde sous pression uniforme 
Nous considérons un tronçon d’une coque parabolique hyperboloïde semi-mince (2a/h=52) 
avec des arêtes droites articulées soumises à une pression normale uniforme. Nous montrons 


les données géométriques et les propriétés mécaniques du matériau sur la figure 3.20. 


LE, 
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À : : L b ; 
L'équation gouvernant la présente coque est donnée par : Z = (2) X.Y. C’est une coque 


distinctive ayant une courbure Gaussienne négative et sa géométrie convient à une coque 
profonde. 
Tableau 3.9. Parabolique hyperboloïde sous pression 


— Caractéristiques géométriques et mécaniques- 


Longueur (2a) 12.92 cm 
Elévation (C) 1.304 cm 
Module de Young (E) 5.10° kg/cm°? 
Coefficient de Poisson (v ) 0.39 
Epaisseur (h) 0.25 cm 


Charge de pression normale (q) -1.0 kg/cm? 


Tous les côtés sont articulés 


La structure n'est pas symétrique. Il est nécessaire de discrétiser l'ensemble de la coque. 
Le but du présent exemple est de voir la capacité de notre élément à modéliser des structures à 


double courbure. Le présent exemple est à flexion légèrement dominante. En effet, au point O 


6Mx 
et pour z=h/2, on a: no 1:5. 


Figure 3.20 Coque Parabolique hyperboloïde sous pression uniforme 
Nos résultats sont présentés sur la figure 3.21. Ils caractérisent la déflexion au centre de 
la coque en fonction du nombre total de degrés de liberté pour différents maiïllages. Nous 
comparons nos résultats avec la solution analytique donnée par Chetty et Tottenham (1964) 
(Déplacement W au centre O égale à 0.00874) et avec les résultats numériques de Guenfoud 
(1990, 1996). Nous notons une bonne convergence monotone des résultats de notre élément 


vers la solution analytique. 
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Hyperbolic paraboloid shell under uniform pressure 
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Figure 3.21 Coque Parabolique hyperboloïde sous pression uniforme 


Flèche au centre de l’hyperboloïde 


3.6.8 Coque hémisphérique (test de Mac Neal) 

Nous étudions la coque hémisphérique soumise à quatre forces ponctuelles radiales en 
auto-équilibre distante à des intervalles de 90 °, deux forces entrantes (vers l'intérieur) et deux 
autres sortantes (vers l'extérieur) aux quarts de ses bords ouverts, via le modèle d’un quart de 
la structure représenté sur la Figure 3.22. Ce problème a pour but de vérifier la performance de 
l'élément à modéliser les rotations de solide rigide et la flexion inextensible pour une coque à 
double courbures. Sur la figure 3.22 nous présentons les propriétés géométriques et matérielles 
de la coque en question. Flugge en (1960) présente une solution analytique au problème. La 
solution de référence pour le déplacement radial aux points chargés est de 0,0924. En raison de 
la symétrie, nous analysons exclusivement le quart de cette coque hémisphérique. Nous 
exposons les résultats pour différents maillages à la figure 3.22. Ils indiquent que l'élément 
proposé fonctionne bien en comparaison avec d'autres éléments de la littérature. Nous 
comparons nos résultats à la solution analytique donnée par Flugge (1960), à la solution de 


référence due à Belytchko (1984) et aux résultats numériques de Guenfoud (1990, 1996). 


Tableau 3.10 Coque hémisphérique (test de Mac Neal) 
— Caractéristiques géométriques et mécaniques- 
“Rayon delasphère(R)  1000m 
Module de Young (E) 6.825x10/Pa 

Coefficient de Poisson (v) 0.3 
Epaisseur (h) 0.04m 
Charge (Concentrée) P=2N 
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Figure 3.22 Coque hémisphérique (test de Mac Neal) 


Données géométrique et maillage 
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Figure 3.23 Coque hémisphérique (test de Mac Neal) 


Déplacement radial au point de chargement 


3.7. Conclusion 

La formulation d'un élément, de coque mince triangulaire à facette plane avec une rotation 
vraie basée, sur l'approche en déformation a été développée avec succès et présentée dans le 
présent chapitre. Les trois déplacements translationnels (u, v et w) sont décrits en termes de 
fonctions polynomiales cubiques. L'utilisation de champs d'ordre égal pour tous les 
déplacements (membranaire et flexionnel) a pour effet de rapprocher plus strictement la 
condition de mouvement du solide rigide. L'originalité dans la formulation du présent élément 
réside dans l'utilisation du modèle en déformation et l'utilisation de concepts et techniques 


suivants pour : 
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- L'enrichissement des champs de déplacement (raffinement P), donc une meilleure 
précision dans l'approximation de la solution. 

- Amélioration du comportement en cas de distorsion géométrique du maillage. 

- L’évitement du problème de blocage de membrane pour les structures courbes. 

- La réponse aux problèmes numériques induits par l'absence de rigidité liée à la rotation 
autour de la normale dans le cas des éléments coplanaires. 

Ces concepts et techniques sont : 

- Adoption de l'approche en déformation ; 

- L'introduction d'un quatrième nœud interne (nœuds esclave) dans l'élément triangulaire 
à trois nœuds ; 

- La réduction des matrices de raideur élémentaires en utilisant la technique de 
condensation statique ; 

- L'utilisation de l'intégration analytique pour évaluer la matrice de rigidité. 

Cette approche nous a offert à un élément de coque à facettes planes, compétitif, robuste 
et efficace. La formulation actuelle (approche en déformation) montre qu'elle est cohérente dans 
une très grande variété de situations d'analyse linéaire. Une série de problèmes test guident à 
évaluer l'efficacité de l'élément par rapport à d'autres éléments se trouvant dans la littérature. 
Les résultats obtenus ont confirmé le taux de convergence rapide de l'élément. L'élément 
proposé a l'avantage d'être de forme simple et utilise les six degrés de liberté. En outre, 1l peut 
être utilisé pour l'analyse de structures à coque mince, même celles avec des géométries 


complexes. 
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ANNEXES 
ANNEXE I - Matrice des coordonnées nodales [A] Elément T43_Eq 
x? : 2 3x2y1 
1101. x! 0 |yil -2vxiyi - 3vxz2 | 6x1y1 AN - 
we y2C24) | -3y2 7 | yo) 
ya s 3Y 2x1 É 2 
0O[1/x 0 !|vilx - -2VX1Y1 | 6xX1Y1 | 3VY 2 - si 
x12(2#v) -3x12 | x1°(2+v) | 
3(yr- AE 
0 0 1 0 0 0 -2x1 2yi 6Y1 -6x1 xp) xp) 
x? 5 ; 3x2y2 
110! |x| 0 | y2 | -2vxy - 3vx22 | 6x2y2 or - 
y? y2(24+v) | -3y2 7 | y#Q+v) 
ya? - ANS |||e: ; 
011] x) 0 | y2)x2 - -2Vx2Y2 | 6x2ÿ2 | 3vy? : Dr 
e x2(2+V) -3x22 | x2*(2+v) 2 
ne 3(y2- | 3(y2- 
0 011010 0 -2X2 2y2 6Y2 -6X2 x) x) (3.53) 
x? à 5 3x2y3 
110! |x3) 0 |y3| -2vxays - 3vx32 | 6x3y3 nent - 
Y3 y#(2+v) | -3y? 3 | y#(2+v) 
ya - 3Y32x3 È 5 
01x30 !|y3|x3 - -2VX3Y3 | 6xX3Y3 | 3vy# - ue 
x32(2+V) -3x32 | x: (2+v) : 
3(y32- 3(y32- 
ol011101010) -2%x: 2: 6y: | -6x: LD D 
x@ E— 3x# ya 
110! [|xal O | yal -2vxaya : 3vx22 | 6x4ÿ4 vor : 
v2(2+v) | -3y2 "| y#@rv) 
y - Ne as 
O1] x4 | O0 | ya) xa - -2VX4Ya | 6X4Ya | 3vy# - Re 
x42(2+4) -3x42 | x4(2+4v) . 
3(y#- 3(y#- 
0 0 1 0 0 0 -2X4 2Y4 6Y4 -6X4 x) x) 
ANNEXE II - Matrice [K, ] Elément T43_Eq 
K,]= fflel DIekxay (3.54) 


H,,= [[x'y'axay (3.55) 
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* Forme générale 


K,]= 


CT. OR © 


1 0 1 
. 1+ v) Jr *É 000010 2y -2-v-x 6-x -6-v:y 6-x-y —6:v:x:y |dxdy 


ml 
! 
D 
< 
< 
NN 
< 
| 
EN 
* 


000002 -4.x.(1+v) 4-y-(1+v) 0 0 —6.x.(1+v) -6.y.(1+v) 


| 0 
0 
| 
| 0 
| 0 
| 0 0 2 000100 -2.v:y 2x —6-v:x 6:y —6:v:x-y 6-x-y (3.56) 
( 
| 
| 0 
| 0 
| 
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* Forme développée avant intégration analytique 


000 0 0 0 0 0 0 0 0 0 
000 0 0 0 0 0 0 0 0 0 
000 0 0 0 0 0 0 0 0 0 
0 0 o — L 0 0 2 0 6 0 6 
x y xy 
(G-) G- 
D0 20-27 | 0 2y 0 6 0 6 0 
X XY 
(=) (1-11) 
2 : : (3.57) 
E..000 0 0 —4 x —4 y 0 0 —6 x —6 y 
(1+ v) 
000 0 2y 4x La v+8x(1+ v) (4 vxy+ 18xy(1 + v)) 12xy -12y° L12xÿ + 12% (1+ v)] [12 vx + 12xÿ (1+ v)] 
000 2x 0 4y  (-vxy+ 18xy1+ v)) [añrsp(irv)] _12VX  12xy lis yé 12 x )) [122 y + 12° (14 v)| 
000 0 6x 0 12 xy _12vé 36  -36vxy 36X y 36 vx y 
000 6y 0 0 -12Y 12 xy 36 vxy  36y° _36 vxÿ 36xŸ 
000 0 és 4e lowivétaul Covsebeétuul: Sy seu bvrénisl) ls re (sl 
[000 6xy 0 6 [vx +122 (14 v)] [iody+2y (ia v)] -36v6y 36x7 [3624182214 v)]l [sé y +18y" (14 v)] 


105 


Chapitre 3 Formulation d’un élément de coque mince triangulaire à facette plane 


* Forme développée après intégration analytique 


K,]= 


000 0 0 0 0 0 0 0 0 0 
0 0 0 0 0 0 0 0 0 
00 0 0 0 0 0 0 0 0 
000,0 YO 0 2H 0 6- Ho 0 6-Hi 
(=) (1-2) 
000 nn _ 0 Pr 0 6: Hio 0 6H 0 
(1-1) (1-2) (3.58) 
E Hoo 
000 0 0  —— 4: Ho 4H 0 0 6H 6: Hp 
(1+ v) 
000 0 2-Hoi —4-Hio  [4:Ho+ 8H: (1+ v)] (4: Hi +18-Hn-(1+v)] 12H -12-Ho  [12-H+12-Ho-(1+v)] [-12:vHp+12-Hn-(1+v)] 
0 00 2-Hy 0 4Hy [4 Hi+18-Hn:(1+v)| C4: Ho + 8: Ho (1+v)] -12-v-Ho 12-Hn [-12-v-Hi+12-Hn-(1+v)] [12H + 2Hal1+ v)] 
000 O0  6-Ho 0 12-Hi _12- Ho 36-Ho  -36-Hy 36: Hi 36: vH 
0 00 6-Hu 0 0 -12- Ho 12-Hi -36- Hi 36- Ho =36- v - Hn 36: Hi 
000 0 6Hn —6-Ho [12-Hp+12-H3o-(1+v)] [-12-v-Hi+12-Hi-(1+v)] 36-Hn -36-v-Hn [36-Hn+18-Hi-(1+v)] [-36-v-Hp+18-Hn-(1+v)] 
[000 6-Hu 0 6: Ho |-12-v-Hp+12-Hp-(1+v)] [12-Ha+12-Hos-(1+v)] -36-v-Hn 36H [-36-v-Hn+18-Hn-(1+v)] [36-H+18:Hy-(1+ v)] 
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ANNEXE III — Formulation des champs d’approximation à partir des solutions homogènes 
(Polynômes bi-harmoniques) de la fonction d’Airy 


IIL.a. Champs des déformations 


N°|F=P;(xy)| 0°F | OF | 0°F 1(0F  Ô°F 1(6F  Ô°F 2(1+v)( °F 
A2 A2 Ex — 4 Mac Pa D au 2 — 
ox OY OxOY E\ ©y ox E\ x ©y E Cx0y 
1 1 0 0 0 0 0 0 
> x 0 0 0 0 0 0 
3 y 0 0 0 0 0 0 
; x + y) V 1 à 1 ; : 
halls E 
2 +y.y? l I 
5 Ê ss ) - £ - 0 0 — 0 
2(1-v ) 1-v 1-v E 
6 0 0 a. 0 0 : 
(i+v) (i+v) Ë 
-2.v. 2. 
7 x°.y 2y 0 2x = à > - Ax(1 +vYE 
2 -2.v. 
8 x.y” 0 2x 2y L = - Ay(1 + vVE 
9 x° 6x 0 0 = ” 0 
10 y° 0 6y 0 y = 0 
- 6.V.X. 6x. 
11 X .) 6x.y 0 3x° . 2 -6x"(1+vyE 
6x. - 6.VX. 
|  xy° o |6xy | 3° _ _ _6y (1+VYE 


(3.59) 
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N Para- Déformation Déplacements Rotation 
° | mètres | E.e E.e, E.y, Eu E.v OM oU E.0. 
d’interpol | de Fe 
a-tion : Ÿ = 1 . av = ou 
2 | ôx Ôy 
1 ai 0 0 0 1 0 0 0 0 
2 a2 0 0 0 0 1 0 0 0 
3 a3 0 0 0 -Y x 1 -] 1 
4 a4 Il 0 0 X 0 0 0 0 
5 as 0 1 0 0 y 0 0 0 
6 a6 0 0 2 y X 1 1 0 
7 a7 -2.V.Y 2.y -Ax(1 +v) -2.V.X.y y” -x/(2+v) -2x(2 +y) -2.v.x -2.x 
8 as 2.x -2vx | - 4y(1 +v) <?- y” (2+v) - 2.V.x.y - 2.V.y - 2.Y(2 +v) 2.Y 
9 a9 -6.v.x 6.x 0 -3.v.x° -3y° 6.x.y 6.y -6y 6.y 
2 2 _ à 

: a10 6.y . 6.v.y 0 6.X.y ee 3.v.y _ 3x 6x 6.x 6.x 
_. -6v.xy| 6xy | -6x°(1+v) -3v.x".y-v" | 3xy°-x°(2+v) 3.y°-3.x/(2+v) -3v.x° -3.y° 3.(°-x°) 
: a12 6.xy |-6V.x.y| -6y”(1+v) 3.x°.y-Y (2+v) -3v.xy”-x° -3v.y"-3x° | 3x°-3y/2+v) 3(ÿ°-x°) 


(3.60 
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ANNEXE IV - ELEMENT FINI DE PLAQUE MINCE 


* Forme générale de [K, | avant intégration 
[K 0 ] = 


© © © 
So © © © 
(e} 


Le ù 0 0 0 1 x y xy 0 0 0 0 0 (3.61) 


2 . 0 0 0 0 0 0 O0 1 x xy 0 |dxd 
{1 v°) xy 0 ( \) Y XY y 
F 2 2 
0 1 0 GES 0 0 0 0 0 2x x O0 2y 0 Y 1 


* Forme développée de [K, | avant intégration analytique 


[O 0 0 0 0 0 0 0 0 0 0 0 
0 0 0 0 0 0 0 0 0 0 0 
0 0 0 0 0 0 0 0 0 0 
| < y X.Y v V.X V.y V.X.Y 0 
x? X.y x2.y v.x vx2 V.X.Y vx2.y 0 
[K | : Eh? à 4x2 + y2 2x + x.y2 V.Y  X.y.(v +44) v.y2 xy2.(v +24) 2xd 3.62 
0 12.(1- v2) x + à + x2.y2 VX.y x2.y.(v +24) x.y? x2.y2.(v+d) x2 à ( 7 ) 
1 X y X.y 0 
Sym. x? + 4y2d X.Y x2.y + 2.3? d 2.y.d 
y? x.y? 0 
x? y2 + y+d y2d 
Æ d | 
* Forme développée de [K, | après intégration analytique 
3 
[K |- nn" 
01 y24-v2) 
[0 0 0 0 0 0 0 0 0 0 0 0 | 
0 0 0 0 0 0 0 0 0 0 0 
0 O0 0 0 0 0 0 0 0 0 
H v.H v.H v.H v.H 0 
00 10 01 11 00 10 01 11 
H H H v.H v.H v.H v.H 0 
20 11 21 10 20 11 21 (3.63 
4H, d+H. 2H, d+H, vH H (v+4d) vH H {(v+2d) 2H, d . 
20 02 30 12 01 11 02 12 10 
H, d+H v.H H,.(v+2d) vH H,.{v+d)  H, d 
40 22 11 21 12 22 20 
H H H H 0 
00 10 01 11 
Sym. H, +4H d H H,,+2H d 2H .d 
20 02 11 21 03 OI 
H H 0 
02 12 
H,,+H d H d 
22 04 02 
d 
L 00 


Avec : d = (=) et Hag = JJ x*yFax. dy 


Chapitre 3 Formulation d'un élément de coque mince triangulaire à facette plane 


ANNEXE V - Intégration analytique sur un triangle 
L'évaluation par intégration numérique (cas des éléments isoparamètrique) sur un 
élément triangulaire (de type Hammer dans ce cas) des différents vecteurs et matrices 
nécessaires pour une analyse quelconque par éléments finis génère les problèmes suivants : 

- Altération de la précision pour les maillages distordus géométriquement, 
particulièrement lorsque les champs d’approximation sont (les polynômes d’interpolations) de 
degré élevé. 

- Le blocage de la solution dans le cas des structures minces. 

L'utilisation de l’intégration analytique était très délicate et pratiquement impossible jusqu’à 
un certain temps, avant l’apparition des machines puissantes et d’outil de calcul utilisant le 
calcul symbolique. Actuellement, l’intégration analytique sur un domaine bidimensionnel 
(triangle et quadrangle) devient possible suite aux logiciels dédiés pour ce but (Mathematica). 


L'évaluation de la matrice [K, |= [Qx») lIDI[Q(xy |dS analytiquement est maintenant 


possible et elle permet d’éviter les problèmes cités en haut. Une démarche simple optée par 
plusieurs auteurs : Onate (2009), Hamadi et al. (2006), Belarbi (2000), Himeur et al. (2008) et 
Guenfoud (2018), consiste en la formulation d’une formule universelle aisé à mettre en œuvre 
et permettant une estimation automatique de la matrice[K, |, ceci quel que soit le degré du 
polynôme résultat du produit [Q,1]=[Qxy/I21[Qx»] et quel que soit la distorsion de 
l’élément. Cette méthode se présente comme suit: Soit un élément fini (triangulaire par 


exemple) (figure 3.24). 


Y 


Ÿ2 
Y3 


Yi 


X1 X3 X2 x 


Figure 3.24 Intégration analytique, Exemple d’un élément fini triangulaire 
Les termes de la matrice résultat du produit [Q,1=[Q(x.y]"12][Q(x.y)] sont sous l’allure 


Cx*.y”. Donc l’intégration de la formule [IQx») l'ID]{Q (x.y ldS revient à intégrer chacun des 


S 


termes de la matrice [Q,], qui se mis sous la forme :- 
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1= [[Cx".y’ ax dy (3.64) 


Avec :I1=11+1 et C : une constante. 
- li: représente l’intégrale entre les droites Y1 et Y3 dans le sens des y et entre x1 et x3 
dans le sens des x. 
- DL représente l’intégrale entre les droites Y et Y2 dans le sens des y et entre x3 et x2 
dans le sens des x. 


Les valeurs de I. et L sont données comme suit : 


x} A 
l;= [fCx".v/ax.d) = [FICx".yPax. dy (3.65) 
A À sh 


Y : représente l’ordonnée de la droite d’équation y = a.x +b. 
donc pour f = 2, nous aurons y? = (a.x + b}? = 1.a2.x2 + 2.a.b.x + 1.b?, 
pour f = 3, nous aurons y° = (a.x + b)° = 1.a°.x* + 3.a2.b.x2 + 3.a.b2.x + 1.b°, 


Ainsi, selon la valeur de B, on aboutit à la forme générale suivante : 
k=B+1 k=B+1 
etre Mit beat = Oba Abr (3.66) 
k=1 KI 
C(k) sont des coefficients dont le nombre et la valeur sont fonction de la valeur de f (voir 
tableau 3.11). 


La sommation continue [y'ay prend la forme : 


[yfay = = y = = | (ax +b}}* = _ | À cutter ait (3.67) 
Tableau 3.11 : Valeur des coefficients C(k) en fonction de f 
CCR), k=1,6 
CD | CO | CG | C@ | CG | CG) 

B =0 il - - - - - 
B=1 il il - - - - 
B =2 il à Il - - - 
B =3 il 3 3 Il - - 
B 4 il 4 6 4 il - 
B =5 il 5 10 10 5 Il 


Pour des bornes d’intégration formées des droites d’équations Yi = ai.x + bi et Yi = a;.x + b;, 


l’équation (3.67) devient : 
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y. 


j k=p+2 
[y'a = Fa: S ONCE ap PTE) ji (3.68) 
Y. k=1 


L’évaluation de l’intégrale de l’équation (3.68) par rapport aux limites d’intégration Xm et Xn 


nous permet d’écrire : 


x (Y x : 
n à Î | =p+2 L : : 3 _ 
[x orale fe > C(k).(a ;* sd ë =" ed Eye l'dx (3.69) 
Xn r; Xm k=1 
L’expression définitive de l’intégrale (3.64) prend la forme : 
1 1 k-1 ET TR k 
Le Fa D pape AAA CIRE Ve DRE 1 GC dE (3.70) 
k=1 


Pour le cas d’un triangle nous avons p=2 


I=Y1, (3.71) 


ANNEXE VI - Condensation statique 


La condensation statique est une technique qui opine de libérer des degrés de liberté 
d’un élément en opérant sur le système d’équation régissant l’équilibre de l’élément (Sarf 
1991), par la simplification de ses composantes (matrice de rigidité, vecteurs de charges). 
L'idée essentielle est d’exprimer les variables inconnues relatives aux degrés de liberté à écarter 
du système d’équations en fonction des variables inconnues à retenir, d’où la réduction du 
nombre d’équations à celui des inconnues restantes. Le privilégie étant la réduction du temps 
de calcul pour la résolution du système d’équations, l’utilisation limitée de la capacité mémoire 
et l’élimination des DDL liés aux nœuds esclaves. Cette technique est une méthode de calcul 
agissant au niveau élémentaire (avant assemblage). Les étapes de calcul sont les suivants. Dans 
la méthode des éléments finis, l’équilibre statique élémentaire est gouverné par un système 
d’équation de type : 

[KAU}={F} (3.72) 
Tel que : [K] : La matrice de rigidité de l’élément, 
{U } : Les déplacements nodaux, 


{F} : Le vecteur chargement. 


Indiquons par l’indice e les équations à éliminer et par l’indice c les équations à garder, puis 


divisons les composantes du système d’équations comme c’est illustré aux tableaux 3.12 à 3.14. 
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Tableau 3.12 : Partition de la matrice de rigidité 


DDL à Garder | DDL à éliminer 


DDL à garder K.] Fe. 


DL Sininér [KR Re) 


Tableau 3.13 : Partition du vecteur des déplacements 


DDL à garder | U.} 


DDL à éliminer | {U,} 


Tableau 3.14 : Partition du vecteur des charges équivalentes extérieures 


DDL à garder {F.} 


DDL à éliminer | {F, } 


Commençant de cette partition, la relation (3.72) devient: 


[K.] [K. U, F 3 73 
[K.] [K. Ü. ù F | 
A partir du système (3.73), nous pouvons soustraire le vecteur des déplacements {U ; } (des DDL 


à éliminer) : 


{U.}= {K.JK.F{U.}-+-[K.F{)} (3.74) 
{0.}=-{K.JIK.F'{U}+ KT} (3.75) 


Le but est d’éliminer les DDL associés au nœud interne. En imposant que la composante du 


vecteur des charges extérieures {F.}=0 et en remplaçant {U,} de l’équation (3.75) par sa 


valeur donnée par l’équation (3.74) nous déduisons : 


{.}= [KR PK K UK Tr} (3.75) 


Cette dernière relation (3.75) se réduit pour livrer une forme semblable aux équations de 


l’équilibre statique dans la méthode des éléments finis : 


UK] K.UK. TK. IKU.} =) (3.76) 
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La matrice [K]= IK..]- K.IK.F K.] forme la matrice de rigidité condensée et le système 
d’équation (3.72) se réduit au système suivant : 
IR NU. }={r) (3.77) 


La relation (3.77) ne contient donc que les paramètres (inconnus) liés aux DDL à garder. 
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4.1. Introduction 


Dans le chapitre précédent nous avons pu développer un élément spécial de coque mince 
triangulaire, à facette plane, à trois nœuds et à six degrés de liberté par nœuds en se basant sur 
la formulation en déformation. La finalité du présent chapitre est d’adapter l’élément développé 
pour qu’il puisse être applicable aux structures coque en matériaux à gradient de propriété 
fonctionnel (FGM). La particularité des structures coque en FGM en comparant avec celles en 
matériaux homogènes isotropes est l’existence d’une matrice de couplage entre le mode de 
déformation membranaire et le mode flexionnel ; ainsi l’intégration sur la coordonnée z n’est 
plus constante suite à la dépendance des caractéristiques mécaniques (E(z) et v(z)) de la 
coordonnée z. Les matrices de rigidités liant les contraintes aux déformations deviennent 


fonction polynomiale de la coordonnée z. 


Toutefois le problème de compatibilité reste poser suite au travail de l’élément membranaire 
de classe C° en contre partie à l’élément flexionnel de classe C! dans la formulation classique 


des éléments finis, ce qui soulève les problématiques suivantes : 


- La pauvreté en déplacement de l’élément membranaire lors de son utilisation à des 
problèmes à flexion dominante, ce qui impose l’utilisation de maillage très fin dans le 
but de se rapprocher de la géométrie curviligne d’une part et de se converger vers la 
solution réelle par un enrichissement de la solution en mode h. 

- La continuité inter élément et la conformité lors de l’utilisation de l’élément 
membranaire de classe C° en conjonction avec l’élément de plaque mince de classe C!. 

- L'absence de la rigidité propre au sixième degré de liberté dans le cas des éléments 


coplanaires. 


L'objectif du présent chapitre consiste à trouver des solutions à ces problèmes où de les 
éviter plus précisément dans le cas de la modélisation des coques en composite (FGM). Nous 
atteindrons cet objectif par le développement d’un élément de coque spécial, pour la 
modélisation des structures en FGM, à facette plane toujours basé sur la formulation en 
déformation. Il est construit par la superposition d’un élément de membrane et d’un élément de 
plaque en flexion tous deux dédiés pour les FGM. A cet effet, une matrice de couplage 
membrane-flexion est nécessaire. Dans ce qui suit nous exposons le développement des 


différentes matrices (membrane, flexion et couplage membrane-flexion). 


117 


Chapitre 4 Elément de coque pour les FGM 


Z 
Céramique 


Figure 4.1 Structure d’une plaque en FGM 


4.2. Formulation de l’élément de membrane pour les FGM 


4.2.1. Elément membranaire «T43_ Eq FGM» 


L’élément visé n’est autre que l’élément de membrane T43 54 développé par Himeur 
(2008) et Guenfoud H. (2018) au chapitre trois dans la section 3.2 généralisé aux matériaux 
composites de type FGM. II est donc triangulaire à quatre nœuds et à trois degrés de liberté par 
nœud (les deux translations u et v et la rotation autour de la normale 6, (Drilling Rotation)) 


dans son repère local, figure (4.2.). 


VA 


Figure 4.2. Elément triangulaire membranaire à quatre nœuds “T43_Eq_FGM” 


Et à trois DDL par nœud (deux translations u et v et la rotation @,) 


Les ingrédients nécessaires pour la construction de la matrice de rigidité du présent 


élément sont : 


- La matrice des coordonnées nodales élémentaire [A] (voir équation (2.15) chapitre 3), 
qu’on va la noter dans le présent chapitre par: [Aras ra rem. Elle dépend 
nécessairement des coordonnées locales x et y des nœuds, du module d’élasticité E(z) 


et du coefficient de Poisson v(z). 


f(x1, V1, z) 

_ 1 f(x, 2,2) Se, 
[Ara3 9 rc] En E(z) F(x3, V3, z) 
f(x4 Ya, z) 


Tel que 
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1 0 -y x 0 y -2vxy a?-y204v) (x + y?) 6xy (x y-y?) Gx2y-y 241) 
[fx v,2)]=10 1 x 0 y x y2-x20+v) -2Vxy 6xy H(vy2 4x2) Gyx-x (24) (y2x-x)) (4.2) 
00 1 000 -2x 2y 6y -6x 42 -x2) 4y2-x?) 


exprimée au quatre nœuds de l’élément. Et : 


A ee (£ F ) LE, 


- La matrice de déformation [Q] donnée par l’équation (3.20) du chapitre 3. Cette matrice 
dépend essentiellement des coordonnées locales x et y des nœuds de l’élément et des 
caractéristiques mécaniques (E(z) et v(z)) du matériau (FGM dans notre cas) constituant 


l’élément. On va la désignée par [Qra3 ra rc |. 


0 0 O0 1 O0 O -2w 2x -6vv  6y - Gvvx 6xy 
Gide 0 0 0 0 1 0 2y 2 6x -6w  6xy _Gwvx 
0 0 0 0 0 2 -4x(l+v) -4y(1+w 0 O0 -6x(1+v) -6ÿ(1+v) 


- La matrice des fonctions de forme exprimé en fonction de la matrice des coordonnées 
[A], de la matrice polynomiale [f(x, y, z)] et du module d’élasticité E(z) du composite 


FGM. Elle est donnée par l’expression : 


1 _ 
[Nras sa ro] = = If, 2)][Aras ra rem] (43) 


E(z) 
Ainsi, pour le cas plan des contraintes et pour les FGM les contraintes sont liées aux 


déformations infinitésimales par la relation suivante : 


1 v 0 
O: € 
. a E(z) v 1 0 (4.4) 
Y(T 1 — v2 1-7 Y 
Try 0 0 ——|Uxy 


Qu'on peut l’écrire sous la forme condensée suivante : 


{o-r43 ka FcM} = [Cras ra rem|{Era3 Eq ru} GS) 


Il est à noter, que dans la présente étude le coefficient de Poisson est pris constant 
(hypothèse généralement prise par plusieurs chercheurs dont Delale (1983), Guenfoud H. 
(2016, 2017) et Ziou (2016, 2017)). 

42.2. Matrice de rigidité membranaire [K%43 a rem] 


Le principe des travaux virtuels pour un élément fini après discrétisation s’écrit : 
P P P P 


(EWint)® = (EWext)° 1:9) 
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Tel que (6W;,)° est le travail virtuel élémentaire interne dépendant des paramètres internes 
(Contraintes et déformations) et (6W,,+)° est le travail virtuel élémentaire externe dépendant 


des paramètres externes (déplacements et forces). 


Le principe s’énonce comme suit : Le travail virtuel interne est égal au travail virtuel 


externe pour tout déplacement virtuel cinématiquement admissible. 


Pour le cas de notre élément, le travail virtuel interne élémentaire s’exprime par une 
intégrale sur le volume élémentaire du produit des déformations virtuelles par les contraintes, 


d’où : 


(SWint)° = [ Den seule [ous ea rem]dV° (4.7) 
ve 
Connaissant que : 
{er43 ta FcM} _ [Nra3 69 rem]{ das ea rem} 


- (4.8) 
{8r43 eq ru} F CE y, 2)| [Arras ea Feu | dise ul 


Et que : 

{ora3 ra FcM} — [Cras ra rom Eras ra rcM} (4.9) 
Et en remplaçant dans l’expression (4.7) {Eras ra rom} et {ora3 ra Fc} par leurs valeurs 
données, respectivement dans les équations (4.8) et (4.9), on obtient : 


(8Wint)° = dd sue | [ [Qra3 ea rem y, NÉ 
ve 


(4.10) 
[Cras ea rc] [Qrassaes y, 2)| [Ass rareml dV°) {dT43 ea FM} 

La matrice de rigidité de l’élément membranaire en FGM s’obtient à partir de 
l’expression entre parenthèse de l’équation (4.10). L'intégration sur le volume de l’élément se 
décompose en une intégration bidimensionnelle sur la surface de l’élément et une intégration 
unidimensionnelle sur l’épaisseur de l’élément. Donc, dV°= dA°.dz=dz.(dx.dy). 

En conséquence, la matrice de rigidité élémentaire pour le cas d’un matériau FGM tirée 
de l’expression (4.10) devient : 


[Résssseul = | [ [Qra3 Ea FeM(X, Y, nes sul 
ve 


(4.11) 
[Cras ea rem] [Qrassares (y, 2)| [ras ea Fe | az. (dx. dy)) 
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4.3. Formulation de l’élément de plaque en flexion pour les FGM 


4.6.1. Elément de plaque mince «PMr43_FGM> 


Au chapitre précédent, nous avons décrit les étapes de construction de l’élément de plaque 
mince (hypothèse de Kirchhoff) triangulaire à quatre nœuds développé par Himeur (2014) et 
Guenfoud (2018). Cet élément est destiné pour un matériau isotrope. La généralisation de 
l’élément pour les FGM est le but du présent paragraphe. L’élément est de forme triangulaire à 
quatre nœuds et à trois de degré de liberté par nœud (deux rotations et une translation 


transversal), figure 4.3. Il est basé sur la formulation en déformation. 


y 


+ Un quatrième nœud fictif 


Figure 4.3 Elément triangulaire de plaque mince 
à trois degrés de liberté par nœud (w, 6, 9,) pour les FGM ( PMr43_FGM) 

Pour la construction de la matrice de rigidité du présent élément de plaque en flexion dédié 
pour les matériaux composites à distribution graduel de ses caractéristiques, nous avons besoin 
des apports suivants ; 

- La matrice des coordonnées nodales élémentaire [A] (voir équation (3.21) chapitre 3). 
On la désigne dans le présent chapitre par : [Apmras reul. Elle dépend nécessairement des 
coordonnées locales x et y des quatre nœuds de l’élément et de la fonction polynomiale 


d’interpolation du champ de déplacement f{(x,y) comme suit : 


f(x1, Y1) 
F2, Y2) (4.12) 
F(X3, Y3) 
FC Y4) 


[Apmras Fc] — 


Tel que 
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L x x” x. X°.y X°.Y y” X.y° y X.Y° X.y 
ND ie Te de ei AT DIT. ct à 
x° x°. ; ‘ 
[fx »]= 0 1 0 x a X.Y 2. 0 . 0 2 : (4.13) 
FO SR Er ee EE A Se ee, Sr € 
2 6 2 2 2 


exprimée au quatre nœuds de l’élément. 


- La matrice de déformation [Q] donnée par l’équation (3.25) du chapitre 3 se dérive 


principalement des coordonnées locales x et y des nœuds de l’élément. Elle est donnée par : 


0 0 O1 x y xy 0 0 0 0 0 
[Qx »]=[0 0 0 0 0 0 0 1 x y xy 0 (4.14) 
0 0 0 0 0 2x x° 0 2y O0 y‘ 1 
Elle est désignée par [ pMTAa3 FGM]. 


- La matrice des fonctions de forme formulée en fonction de la matrice des coordonnées 
nodales [A] et de la matrice polynomiale [f(x, y)]. Elle est donnée par l’expression (3.23) du 
chapitre 3 : 

-1 
[Npuras rem] = [f(& Y)]I[Apmras rom] (4.15) 
Ainsi, pour le cas plan des contraintes et pour les FGM les contraintes sont liées aux 


déformations infinitésimales par la relation suivante : 


d2w 
Ox 1 v O0 Ex 1 v O0 _. 
(e] 
Cy _ v 1 A de ce) v 1 s A (4.16) 
t FO 20 A 0 
XY 2 y 2 d2w 
ôx0y 


P 
Qui dépend de la coordonnée z à travers :  E(z) = (Ec — En) É + :) + En. 


En terme de relation "Moments - Courbures" cette relation peut s’écrire sous la forme 


suivante : 
— 92w 
ox? 
Et Sr 
My = D 1 y Ky — D 1 — V dy? (4.17) 
Myy 0 0 g Kxy 0 0 g re 
Ôx 0y 
Avec 
1 +h/2 
D = —— E(z)z’dz (4.18) 
ES) —h/2 
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Il est à noter, qu’ici dans notre étude le coefficient de Poisson est pris constant 
(hypothèse généralement prise par plusieurs chercheurs dont Guenfoud (2016, 2017) et Ziou 


(2016, 2017))et  E(z) = (E-—En) É £ à +E. 


4.3.2 Matrice de rigidité élémentaire [Kéras rcm| 


Le travail virtuel élémentaire intérieur discrétisé s’exprime comme suit : 
T 
(8Wint)° = [ 6{EpmTas FGM} {OpmTa3 FcM}AV* (4.19) 
ve 


Nous savons que : 


{EpmTas rcM} = [NomrTas reM]{dSmras FcM} 
(4.20) 


—1 
{EpmTa3 FGM} — [QpmT43rçu y, 2)][Apmras Fc] {dE mra3 Fc} 
Et 
{opmras.FeM} = [Cpmra3 FM ]{EpmT43 FcM} (4.21) 
Par substitution dans l’expression (4.19), (Epmras rc m} êt {Opmras Fc m} par leurs 
formules données en (4.20) et (4.21) respectivement, nous obtenons : 
e e T T 1 T 
(OWint) Æ Of Dnra3reu) [OpmTa3rem ( y)] [ApmTasre ] [Comrasrew] 
V (4.22) 
1 
[Qpmras rem (x, y )| [Apmras re] ave) (Qémras rcu} 
De ce fait, la matrice de rigidité élémentaire est représentée par l’expression entre parenthèse 
de l’équation (4.22), elle s’écrit comme suit : 
2 T 7 
[Kémras Feu] = [QpmTa3re (x, y)| [A pmTa3r em ] [Cymrasreu] 
V (4.23) 
1 
[Qomras rem (x, y )] [Apmras Fou] ave) 


1 v O0 
ju E(z)z?dz 
— =} 
Avec el ee É : Se 


2 


L'évaluation de cette matrice de rigidité est établie par intégration analytique des 


différentes composantes, dont les expressions sont de la forme H,,8 = C xtyb. 


Finalement la matrice de rigidité élémentaire de la plaque mince FGM à considérer à 


l’assemblage des différentes matrice pour construire la matrice de rigidité globale de la 
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structure, est celle obtenue après condensation de la matrice Rae ï] pour éliminer les 


degrés de libertés du quatrième nœuds (nœud esclave). 


4.4. Formulation de l’élément de couplage membrane-flexion pour les FGM 


Au chapitre trois dédié au matériau isotrope n’apparait pas une matrice de couplage entre 
le mode de déformation membranaire et celui de flexion suite aux raisons suivantes : 

- Le matériau est isotrope, donc le plan moyen de la coque se trouve à m1- portée entre 
les deux faces extrêmes de l’élément de coque (couplage dû au comportement 
mécanique). 

- L'élément fini de la coque est supposé à facette plane, ce qui élimine le phénomène de 
couplage entre la membrane et la flexion pouvant apparaitre si l’élément est courbe 
(couplage dû à la géométrie courbe). 

Dans le cas de la modélisation des coques en FGM, bien que l’élément considéré soit à 
facette plane, un couplage entre les deux modes de déformations (membrane et flexion) se 
produit en créant ainsi une matrice de couplage membrane flexion qui n’apparut pas dans le cas 
des matériaux isotropes. La raison de ce phénomène revient à l’existence de la coordonnée z 


dans les expressions des déformations (membrane et flexion) avec différents exposants (paire 


D nr : —h 
et impaire) d’une part et les bornes d’intégrations sur z qui ne sont plus entre + = 


Rigidité de couplage membrane-flexion 


La rigidité de couplage membrane-flexion vient du travail virtuel de l’effet de la 
membrane sur la flexion d’une part et de l’effet de la flexion sur la membrane d’autre part. En 
effet, ce travail est le résultat du produit de la déformation virtuelle de membrane multiplié par 
le tenseur de contrainte due à la flexion en plus du produit de la déformation virtuelle de la 
flexion multiplié par le tenseur de contrainte due à la membrane. 


Ce qui peut être traduit par : 


(SWint)° = [ beau [opmras rem]dV® 
ve 


(4.24) 
T 
+ [ {EpmTa3 FcM} [ora3 59 rem]dV® 
ve 
AVE : 
{6843 ra rem} = [Nras ea rem|{8dT43 ea Feu) 
L (4.25) 
{5ëras 89 rem} = [Qra3 £a reM]lATa3 £a rem] {843 ka rc} 
Et 
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{OpmTa3 rcM} = [Cpmras Fam Epmra3 FcM} (4.26) 
Ainsi, 
{8EpmTas FGM} L [NomTas FGM Sd mras rem} 
2j (4.27) 
{8Epmra3 FM} L [Qpmras FGM]|ApmTa3 FGM] {5 mras FM} 
Et 
{Ora3 ra rom} = [Cra3 ra rom {Era ra rom} (4.28) 


En substituant les expressions (4.25) — (4.28) dans la relation donnant le travail virtuel 


interne élémentaire, nous obtenons ce qui suit : 


. T 
(SWint)° = [ ([Qras 84 romÎlAras ra rem] {6043 ra rem} [Cpmra3 rc] 
ve 


-1 
([Qomra3ren) [ApmTa3rçu] {dSmrasren)) dv° 


(4.29) 
1 T 
Ÿ je ([Qpmrsron] [ApmTa3ren | {SaSmrasreu)) [Craszqeu] 


([Qrss ra rom] [Aras ea Feu | {das pa reu}) dv 


L'intégrale sur le volume se décompose en une intégrale sur la surface du triangle dA et une 


intégrale unidimensionnelle sur z. Nous tirons l’intention du lecteur que : 


[ Cras Eq FGM : Z- dz (4.30) 
Z 
Et 
[ Ion ro . Z dz (4.31) 
Z 
Deviennent : 
f : 1 v 0 
E(z)zdz 
[ RUE. (4.32) 
LT =:v2 1—-v 
Z 0 0 g 


Qui représentent la rigidité des couplages membrane-flexion et flexion-membrane et 


E(z) = (EE) É F. à] Es: 


La matrice de couplage : 


RE el 


résultante de l’équation (4.24) fait intervenir les degrés de liberté de membrane et de flexion. 
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4.5. Formulation de l’élément de coque pour les FGM 


Une fois les trois matrices de rigidité (membrane, flexion et couplage) sont évaluées, il 
nous reste qu’à combiner ces trois matrices pour formuler la matrice de rigidité de l’élément 
de coque triangulaire à facette plane. 

Les ingrédients utilisés pour arriver à ce but sont : 

- Utilisation de l’élément triangulaire T43_Eq_FGM pour ce qui est de la rigidité de membrane. 
- Utilisation de l’élément triangulaire de plaque mince PMT43_FGM pour ce qui est de la 
rigidité de flexion 


- Utilisation de l’élément triangulaire de couplage membrane-flexion (T43cmf_FGM). 


La combinaison des trois matrices permet de conduire à la construction de la matrice de 
rigidité de l’élément de coque triangulaire pour les FGM. Cette matrice est pleine et symétrique. 
L’élément de coque peut avoir une orientation quelconque dans l’espace et dans le repère global 


d’une structure XYZ ; et il se représente dans son repère local xyz comme sur la figure (4.4). 


Cal 


Figure 4.4. Représentation de l’élément de coque pour les FGM par rapport aux deux 
repères (global et local) 

L'évaluation de la matrice de rigidité de l’élément de coque pour les FGM dans le repère 

global d'assemblage s’accomplit au biais de la matrice de rotation [RO] comme suit : 
(OU)IKEKU} = (U)IRONTIKÉTIROIU} = (ôu)[KFI{u} 433 

[K6] = [ROFTKŸITRO] 

et cela après éclatement des termes de rigidité des trois modes (membrane, flexion et couplage 
membrane-flexion) à l’intérieur de la matrice de rigidité de l’élément de coque de taille (18x18) 


au niveau élémentaire avant d’entamer l’opération d’assemblage qui se schématise sur les 


figures 4.5 et 4.6. 
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Elément de membrane 
3 
DE (u, v, au, v, 8. u, v, 6, JT 
l 2 
Elément de flexion 
{wi Gus Os We Ge Gps Wa One Oye } 
9DDL 
2 
Elément de couplage 
3 
U, Vi Wi Os Op: Où Us Vs We Oxz Ov O2 à 
18DDL U; V3 Wa 0x2 ya 022 
1 2 
Elément de coque 
L 
Us Vi Wi rs Op: Ou Us Vs We Oxz Op O2 
I8DDL U; V3 Wa 0x2 ya 022 
1 2 


Figure 4.5 Construction de l’élément de coque pour les FGM 
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Dai 
Lt. 2 


HN Terme derigidité de couplage 
% Terme dela rigidité membranaire 
Terme dela rigidité flexionnelle 


Figure 4.6 Structure type de la matrice de rigidité de l’élément coque 
FGM dans le repère local 


4.6. Validation de l’élément de coque pour les FGM 


Afin de mettre en valeur l’élément développé vis-à-vis des deux comportements de base, 
à savoir le comportement membranaire et le comportement flexionnel, nous avons traité les cas 
tests en dessous. Le résultat obtenu est comparé, à la solution de référence se trouvant dans la 


littérature. 
4.6.1. Validation vis-à-vis du comportement flexionnel 


Exemple 1. Console en FGM soumise à la flexion pure. 


Afin de voir la validité de notre élément à modéliser correctement le comportement 
flexionnel d’une structure en FGM; nous considérons la console en FGM de la figure 4.7. Elle 
a une épaisseur h=0.1m, une longueur L=10m et une largeur b=1m. Elle est composée d’un 
mélange de deux matériaux différents, à savoir du métal et de la céramique ; et elle est sollicitée 


à la flexion par une charge P=0.1Pa (figure 4.7.). La surface supérieure est faite entièrement en 
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céramique tandis que celle inférieure est en métal. Les propriétés mécanique de la céramique 
2 E 2 
sont Eç = 1.2E6 GPa,v. = 0 et celles du métal sont E,, = = Vm = (0. Les coordonnées x et y 


sont dans le plan de la console tandis que l’axe z est perpendiculaire au plan de la console 


(dirigé selon l’épaisseur). 


y.  d 


Figure 4.7. Une console en FGM soumise à la flexion par une charge ponctuelle 


Tableau 4.1. Les déplacements w au point de coordonnées (10,0.5) en fonction des 
différents paramètres du matériau 


Paramètre du Ziou (2017) Present element Solution de 
Matériau (P) W U W U comparaison 

0 3,42E-01 | 0,00E+00 | 3,40E-01 | -0,00E-04 3,33E-01 

l 4,57È-01 |-3,86E-04| 4,60E-01 | -3,90E-04 4,62E-01 

2 4,89E-01 | -4,68E-04| 4,92E-01 | -4,70E-04 4,99E-01 

3 5,07E-01 | -4,64E-04| 5,11E-01 | -4,67E-04 5,15E-01 

4 5,21E-01 |-4,39E-04| 5,23E-01 | -4,43E-04 5,25E-01 

. 5,33E-01 | -4,09E-04| 5,34E-01 | -4,13E-04 5,34E-01 


Déplacement w au point de coordonnées (10,0.5) en 
fonction des différents paramètres du matériau 


6.00E-01 


5.00E-01 

4.00E-01 

3.00E-01 

2.00E-01 

1.00E-01 

0.00E+00 
0 1 2 3 4 5 


BZiou MPresent 8 Reference 


Figure 4.8. Les déplacements w au point de coordonnées (10,0.5) en fonction des 


différents paramètres du matériau 
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Le tableau 4.1 et la figure 4.8 regroupent les résultats du déplacement w à l’extrémité 
de la console en FGM (x=10,y=0.5) en fonction des différents paramètres du matériau. Nous 
comparons nos résultats avec ceux de Ziou (2017) et la solution de référence (solution 
analytique). Nous constatons que les déplacements de la console complètement en céramique 
sont moins importants que ceux de la console en FGM. Ceci est dû au fait que les valeurs les 
plus élevée de p correspondent à une haute contribution du métal. Nous notons aussi que les 
déplacements transversaux augmentent quand l’indice du matériau p augmente. Les résultats 
obtenus avec le présent élément sont très performants en les comparants avec ceux obtenus par 
Ziou en 2017 et la solution de référence. Nous notons aussi l’existence du déplacement axial 
non nul, dans le cas de la console en FGM, bien que le test soit à flexion pure. Chose qui 


n'apparait pas dans le cas de la console homogène isotrope (p=0). 


Exemple 2. Une plaque carrée en FGM simplement appuyée soumise à une charge 
concentré ou répartie 


Considérons une plaque carrée en FGM composée d’une part de la céramique (Alumine) 
sur la face supérieure et d’autre part du métal (Aluminium) sur la face inférieur. Les 
caractéristiques mécaniques des deux matériaux sont : Céramique (alumine) : E. = 380 GPa, 
métal (aluminium) : E,, = 70 GPa, v. = v,, = 0.3. La plaque est simplement supportée et elle 
est sollicitée une fois par une charge concentrée au milieu F et une autre fois par une charge 


uniformément répartie qo. 


Céramique 


Figure 4.9 Plaque carrée en FGM simplement appuyée soumise à la flexion 

Les résultats des déplacements adimensionnels au centre (w.) de la plaque FGM 
soumise à une charge concentré F pour un premier cas et une charge répartie g, dans un 
deuxième cas se normalisent respectivement par les équations suivantes : 
Me? 


ar X 10 D= oi E(z)z‘dz (4.34) 


w.D 
—— x 10 
Fa2 


W. : Le déplacement au centre de la plaque FGM, 


a : La longueur du côté de la plaque. 
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Les deux tableaux 4.2 et 4.3 et les figures 4.10 et 4.11 rassemblent les déplacements 


adimensionnels au centre de la plaque FGM pour : les deux cas de chargement, les différentes 


valeurs du paramètre du matériau et pour différents maillage. Les résultats obtenus par le 


présent élément sont en bon accord avec ceux données par Shi-Rong Li en 2015 et Ziou en 


2017. Ils démontrent une convergence rapide de notre élément vers les solutions de référence 


pour les deux cas de chargement et pour les différentes valeurs du paramètre du matériau p. 


Tableau 4.2. Les déplacements adimensionnels W. au centre d’une plaque carrée en FGM 


simplement appuyée soumise à une charge concentrée 


Paramètre | Maillage 2x2 4x4 8x8 16x16 Shi-Rong Li 
P 
Ziou 2,0245 | 2,1035 2,1365 2,1392 21371 
P=0 Présent | 2,0327 | 2,1187 2,1367 21972 
P=1 Ziou 3,4527 | 3,7243 3,7822 3,7794 4,2876 
Présent | 3,4599 | 3,8376 3,9831 4,2754 
P=2 Ziou 4,0253 | 4,4402 4,5049 4,4961 5,4947 
Present | 4,1478 | 4,5712 4,8495 5,3137 
P=3 Ziou 4,3770 | 4,8400 4,9097 4,8995 6,0370 
Present | 4,9658 | 5,2643 5,8219 6,0083 
P=10 Ziou 5,9299 | 6,3482 6,4480 6,4459 7,1333 
Present | 6,0376 | 6,7965 6,9031 7,0765 
Métal Ziou 11,0012 | 11,4187 | 11,5985 11,6129 11,602 
Present | 11,0064 | 11,4389 | 11,5988 11,6123 
Déplacement adimensionnel Wc (p=0) Déplacement adimensionnel Wc (p=1) 
16 4,4 
14 a 4,2 
2,12 4 
" 3,8 
2,06 is 
2,04 34 
2,02 3,2 
2 3 
0 5 10 15 0 10 15 20 
——Ziou —e—Present Reference ——Ziou —e—-Present Reference 
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Déplacement adimensionnel We (p=2) Déplacement adimensionnel We (p=3) 
6 6,5 
5,5 6 
à = 55 a 
5 
5 
4,5 A  — = 
4,5 
L 4 
3,5 35 
3 3 
0 5 10 15 20 0 5 10 15 20 
—e—Ziou —e—Present Reference —e—Ziou —e—-Present Reference 
Déplacement adimensionnel Wc (p=10) Déplacement adimensionnel Wc (Métal) 
7,5 11,7 
11,6 — 
( 115 
6,5 11,4 
11,3 
6 11,2 
& 11,1 
5,5 
11 
5 10,9 
0 5 10 15 20 0 5 10 15 20 
—e—Ziou —e-Present Reference —e—Ziou —e-Present Reference 


Figure 4.10. Les déplacements adimensionnels au centre d’une plaque carrée en FGM 
simplement appuyée soumise à une charge concentrée en fonction du paramètre du matériau p. 


Tableau 4.3. Les déplacements adimensionnels W. au centre d’une plaque carrée en FGM 
simplement appuyée soumise à une charge répartie 


Paramètre P | Maillage | 2x2 4x4 8x8 16x16 | Shi-Rong Li 

Ziou | 0,5061 | 0,6985 | 0,7391 | 0,7408 0,7483 
P=0 Présent | 0,5074 | 0,6995 | 0,7398 | 0,7488 

P=1 Ziou | 0,8631 | 1,2316 | 1,2968 | 1,3966 1,5012 
Présent | 1,0876 | 1,3258 | 1,4721 | 1,5008 

P=2 Ziou 1,0063 | 1,4640 | 1,5357 | 1,8701 1,9238 
Present | 1,1954 | 1,4710 | 1,5892 | 1,9106 

P=3 Ziou 1,0942 | 1,5953 | 1,6727 | 1,8976 21137 
Present | 1,2745 | 1,7019 | 1,8792 | 2,067 

P=10 Ziou 1,4824 | 2,1012 | 2,2150 | 2,3187 2,4976 
Present | 1,5801 | 2,1976 | 2,3391 | 2,4892 

Métal Ziou |2,7475 | 3,7920 | 4,0122 | 4,0439 4,0620 
Present | 2,8460 | 3,8723 | 4,0276 | 4,0641 
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Déplacement adimensionnel We (p=0) 
0,8 
0,75 
0,7 
0,65 
0,6 
0,55 
0,5 
0,45 
0,4 

10 1 


e 
un 
un 


—e—Ziou —e—-Present Shi-Rong Li 


Déplacement adimensionnel We (p=2) 


1,8 
1,6 
1,4 
1,2 
1 
0 5 10 15 
—e—7Jiou —e—-Present Shi-Rong Li 
Déplacement adimensionnel We (p=10) 
2,2 
2 
1,8 
1,6 
1,4 
1,2 
1 
0 5 10 15 


——Ziou —e—-Present Shi-Rong Li 
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Déplacement adimensionnel We (p=1) 


0 5 10 15 20 
——Ziou —e—Present Shi-Rong Li 
Déplacement adimensionnel We (p=3) 
0 5 10 15 20 
——7iou —e-Present Shi-Rong Li 
Déplacement adimensionnel We (Métal) 
0 5 10 15 20 


—e—Ziou —e—Present Shi-Rong Li 


Figure 4.11. Les déplacements adimensionnels W. au centre d’une plaque carrée en FGM 
simplement appuyée soumise à une charge répartie en fonction du paramètre du matériau p. 


Nous remarquons aussi que nos résultats et ceux de Ziou (2017) sont proches les uns aux autres 


pour le cas de matériaux isotrope et pour le cas de matériaux FGM avec p paire. Nos résultats 


pour des p impaires sont plus performants que ceux de Ziou (2017). 


4.6.2. Validation vis-à-vis du comportement membranaire 


Nous reprenons l’exemple précédent de la console sous flexion avec les mêmes données 


géométriques et mécaniques du matériau ; mais cette fois c1 la console FGM est sollicitée à la 


membrane pure (traction pure) par une charge P=0.1Pa, comme s’est représentée sur la figure 


4.12. 
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Figure 4.12 Console en FGM soumise à la traction. 


Tableau 4.4 Les déplacements u au point de coordonnées (10,0.5) en fonction des différents 


paramètres du matériau 


Paramètre du Ziou (2017) Present element | Solution de 

Matériau (P) U W U W  |comparaison 
0 8,59E-06 | 0,00E+00 | 8,47E-6| 0,00 8,33E-06 
1 1,15E-05 | 3,82E-04 |1,15E-5 |3,83E-4) 1,15E-05 
2 1,29E-05 | 4,64E-04 | 1,30E-5 |4,65E-4| 1,31E-05 
3 1,37E-05 | 4,60E-04 | 1,38E-5 |4,61E-4| 1,39E-05 
4 1,43E-05 | 4,34E-04 | 1,44E-5 |4,35E-4|  1,44E-05 
5 1,47E-05 | 4,05E-04 | 1,47E-5 |4,19E-4) 1,47E-05 


Déplacement U au point de coordonnées (10,0.5) en 
fonction des différents paramètres du matériau 


2.00E-05 


1.50E-05 


F1 M 
1.00E-05 F 
5.00E-06 | | 
0.00E+00 u | = 
0 1 2 3 4 5 


BZiou MPresent Reference 


Figure 4.13 Les déplacements u au point de coordonnées (10,0.5) en fonction des différents 
paramètres du matériau 

Le tableau 4.4 et la figure 4.13 regroupent les résultats du déplacement U à l’extrémité 

de la console en FGM (x=10,y=0.5) en fonction des différents paramètres du matériau. Nous 


comparons nos résultats avec ceux de Ziou (2017) et la solution de référence (solution 
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analytique). Nous relevons que les déplacements U de la console entièrement en céramique sont 
plus petits que ceux de la console en FGM. La raison c’est que les valeurs les plus élevée de p 
concordent à une forte contribution du métal. Nous notons aussi que les déplacements 
membranaires U accroissent quand l’indice du matériau p augmente. Nous notons pareillement 
l’existence du déplacement transversal W non nul, dans le cas de la console en FGM, bien que 
le test soit à membrane pure. Chose qui n’apparait pas dans le cas de la console isotrope 
(p=0).Les résultats obtenus avec notre élément sont très performants en les comparants avec 


ceux obtenus par Ziou en 2017 et la solution de référence. 


4.7. Conclusion 


L'objectif principal du présent chapitre est la construction d'un élément fini de coque 
FGM permettant de prendre en compte l'essentiel des situations rencontrées dans le calcul des 
ouvrages. Nous avons développé un élément plat à facette plane triangulaire à trois nœuds et à 
six degrés de libertés par nœud, avec un vraie sixième DDL de rotation autour de la normale, 


permettant la discrétisation des structures coques FGM. 


Le développement de notre élément est basé sur une formulation dite spéciale fondée 
sur l’approche en déformation, ce qui a permis l’enrichissement et le haussement du degré de 
polynôme d’interpolation des deux champs de déplacements, membranaire et flexionnel, avec 
le même degré. Cette pensée et ce raisonnement permet d’éviter le problème d’incompatibilité 
entre les déplacements membranaires et flexionnels. 

La bonne performance de notre élément est clairement observable à travers la série de 


tests de validation. 
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* Forme développée avant intégration analytique 
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ANNEXE II - Matrice [K, | Elément de plaque mince FGM ‘PMras rc” 
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[Ko | = Il Re 2 RC 2. ldxdy | dz (441) 
(=v) x'd+x"y" vxy x“y(v+2d) vxy” x‘y‘(v+d) xd 
1 x y XY 0 
x? +4dy?  xy  x2y+2dy*? 2yd 
y° xy° 0 
x2y?+y#d y’d 
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* Forme développée après intégration analytique 
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Hoo Ho Hoi Hi: 
Ho H:1 H2: 
fE(z)z?dz A4H2od + Hoz  2H3od + H23 
[Ko ] : A=v2) H4od + H22 


ANNEXE III - Matrice [K, | Elément de couplage membrane-flexion 
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* Forme générale (Couplage membrane-flexion) 


cpm-f 
[Ke] 
0 0 0 
0 0 0 
0 0 0 
1 0 0 
0 L 0 1 v O0 
_ JEG)z [= : ee. v 1 0 
an] 2 2%  -KG+v 1—v 
( ) 2 2X —2vx —4y(1+v) 0 0 2 
—6VX 6x 0 
6y —6Vvy 0 
—6vxy  6xy —6x?(1+1v) 
6xy  —6vxy —6y?(1+7v) 
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0 0 
0 0 
0 0 
vH:: 0 
VvHb: 0 
Hi(v+24) 2Hyd 
H(v+d)  Hzod LA 
H;: 0 
H:1 + 2dH,3 2Hoid 
H;2 0 
H>> + Houd Hood 
Hood 
0 0 0 (4.43) 
—2y O0 —y?2 —1|dxdy | dz 
Y  Xy 0 
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* Forme générale (Couplage flexion-membrane) 
cpf-m 
Ko” "| 


0 0 0 
0 0 0 
0 0 0 
1 0 0 
ds 0 1 0 
[E@z RER RS 0 0 2 (4.44) 
er Il 0 O0 O0 O0 O0 —2x —x2 O0 —2y O —y2 —1 1=vlEgl 2 2y —4x(l+v) |dxdy | 4 
0 0 0 0 0 0 0 1 x y xy 0/0 0 —— 2X —2vx —4y(1+v) 
—6Vx 6x 0 
6y —6vy 0 


—6vxy  6xy  —6x?(1+1v) 
6xy  —6vxy —6y?(1+7v) 
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1. Conclusions générales 

Le but principal de cette thèse est de représenter une contribution à la modélisation par 
éléments finis des structures (poutre, plaque et coque) en matériau à gradient fonctionnel. 

Au premier volet de cet recherche, nous avons développé une nouvelle fonction 
polynomiale de cisaillement transversal basé sur une théorie de poutre à déformation de 
cisaillement d'ordre supérieur, satisfaisant les conditions aux limites de contraintes de 
cisaillement nulles (conditions aux limites exactes sur la contrainte de cisaillement) à la fois sur 
les surfaces supérieure et inférieure d’une poutre FGM. La fiabilité, la précision, la performance 
et la polyvalence de la fonction de cisaillement développée ont été évaluées à travers une 
conséquence d’applications en flexion des poutres FGM avec différents cas de chargement, de 
géométrie et des conditions aux limites. Nous avons comparés nos résultats avec des solutions 
analytiques de références (Benchmark) et ceux obtenus par d’autres modèles. Pour les 
différentes analyses statiques les exemples numériques montrent que la fonction développée est 
capable de donner d'excellents résultats. 

L'objectif du deuxième où le second volet de la présente recherche est la mise au point 
d'un élément fini, de coque mince à facette plane basé sur la formulation en déformation, 
triangulaire à trois nœuds et à six degré de liberté par nœud, destiné à la modélisation des 
structures courbe en FGM et permettant de prendre en compte l'essentiel des situations 
rencontrées dans le calcul des ouvrages. Nous avons adopté une approche basé sur l’utilisation 
d’élément à facette plane pour discrétiser les structures coques FGM. Le privilégie fondamental 
des éléments appartenant à cette classe est leur simplicité proportionnelle. En effet leur 
formulation repose sur les théories de membrane et des plaques. De ce fait, on contourne toutes 
les complexités dues à la prise en compte des courbures dans les théories des coques curvilignes. 
Par conséquent, la géométrie des coques FGM peut être approchée en utilisant les éléments 
plans par superposition d’un élément de membrane, d’un élément de flexion et un autre élément 
de couplage (dans le cas des FGM). 

La comparaison des résultats obtenus avec des solutions de références, a montré la 
performance et la précision de l’approche proposée. 

Nous notons que le présent travail de développement a été assidu par un autre de 
programmation en vue d’intégrer ces éléments au code de calcul « Eléments finis » disponible 


au niveau du Laboratoire de Génie Civil et d’'Hydraulique de l’université de Guelma. 
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La programmation est faite avec le langage Fortran d’une part et le progiciel 
Mathematica d’autre part, et elle a consisté à intégrer les présents développements au code de 


calcul d'éléments finis du LGCH. 


2. Perspectives pour l’étude future 

La maturation de l’élément de coque en FGM reste nécessaire. En effet, un post- 
processeur de détermination des contraintes en chaque point et à chaque côte doit être mis au 
point. Ainsi, pour élargir l’application du présent élément de coque FGM aux analyses 
dynamique et de stabilité, nous devons développés les matrices de masse et des contraintes 


initiales ; chose qui n’est pas facile. 
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